Potassium limitation of forest productivity – Part 2: CASTANEA-MAESPA-K shows a reduction in photosynthesis rather than a stoichiometric limitation of tissue formation

Potassium (K) availability constrains forest productivity. Brazilian eucalypt plantations are a good example of the K limitation of wood production. Here, we built upon a previously described model (CASTANEA-MAESPA-K) and used it to understand whether the simulated decline in C source under K deficiency was sufficient to explain the K limitation of wood productivity in Brazilian eucalypt plantations. We developed allocation schemes for both C and K and included these in CASTANEA-MAESPA-K. Neither direct limitations of the C-sink activity nor direct modifications of the C allocation by K availability were included in the model. Simulation results show that the model was successful in replicating the observed patterns of wood productivity limitation by K deficiency. Simulations also show that the response of net primary productivity (NPP) is not linear with increasing K fertilisation. Simulated stem carbon use and water use efficiencies decreased with decreasing levels of K availability. Simulating a direct stoichiometric limitation of NPP or wood growth was not necessary to reproduce the observed decline of productivity under K limitation, suggesting that K stoichiometric plasticity could be different to that of N and P. Confirming previous results from the literature, the model simulated an intense recirculation of K in the trees, suggesting that retranslocation processes were essential for tree functioning. Optimal K fertilisation levels calculated by the model were similar to nutritional recommendations currently applied in Brazilian eucalypt plantations, paving the way for validation of the model at a larger scale and of this approach for developing decision-making tools to improve fertilisation practices.

1 Introduction

Fertilisation trials in tropical eucalypt plantations have been conducted over multiple rotations (Laclau et al., 2010; Battie-Laclau et al., 2016; Gazola et al., 2019; Gonçalves, 2000). These experiments have shown that nutrients can strongly affect tree growth in these highly productive stands. These nutrient limitations can be explained in part by low nutrient supplies from highly weathered soils and in part by the large exports of nutrients with trunk wood (every 6–7 years) from the stands. The frequent export of trunk wood in these fast-growing plantations leads to the export of massive amounts of nutrients that are immobilised in wood (Cornut et al., 2021). In commercial plantations, this issue is solved through the use of fertilisers (NPK, dolomitic lime, and micro-nutrients).

Potassium (K) has been identified as the most limiting nutrient for wood productivity in many omission trials (Gonçalves, 2000; Rocha et al., 2019). Potassium nutrition impacts wood growth through different physiological mechanisms that have been reviewed in detail (Sardans and Peñuelas, 2015; Cornut et al., 2021). In brief, K deficiency is known to depress the assimilation of carbon by the plant (C-source processes, gross primary production), as well as the allocation and use of carbon by the plant for growth (C-sink processes, net primary production).

The study of the limitation of gross primary production (GPP, C source) by K deficiency was explored at the stand level in Part 1 (Cornut et al., 2023) of this series of two papers, using a coupled C–H₂O–K mechanistic model (CASTANEA-MAESPA-K). The simulations showed a strong response of GPP to K deficiency. The GPP in the simulated K omission stand was less than half of that in the simulated fertilised stand (Table 2, Cornut et al., 2022). These results were consistent with previous measurements (Epron et al., 2012) and modelling work (Christina et al., 2018). The strong response of GPP to K availability was due to a reduction in leaf area index and in the photosynthetic capacity per unit area of leaf. The reduction in canopy area was due in part to a slight reduction in leaf production, in part to a decrease in individual leaf area, and in part to a strong decrease in leaf lifespan (Fig. 4, Cornut et al., 2022). The reduction in the photosynthetic capacity of the canopy was associated with the appearance of leaf symptoms in the K-deficient (oK) stand. The impact of symptom area on leaf photosynthetic capacity was sufficient to explain most of the reduction in leaf-scale assimilation in the unfertilised case. A decrease of WUE_GPP (ratio of GPP to transpiration) in the simulated oK stand was also simulated. A sensitivity analysis of the model parameters showed that competition between the organs (trunk, branches, bark, and roots) and leaves for K access had an important impact on GPP in the oK stand (the leaf-to-phloem resorption resistance, R_{leaf→phloem}; Fig. 7 in Cornut et al., 2022). This underlined the need for a precise understanding of K circulation and stoichiometry in the plant.

The impact of K deficiency on net primary productivity (NPP) is a widely observed phenomenon in non-planted forests (Tripler et al., 2006; Baribault et al., 2012) and in planted-forest omission experiments (Laclau et al., 2010; Battie-Laclau et al., 2016; Gazola et al., 2019). While a negative effect of K deficiency on GPP is well documented, the magnitude of its contribution to the decrease in trunk NPP remains unclear. Potassium deficiency also impacts phloem mobility (Epron et al., 2016; Marschner et al., 1996), as well as the loading and unloading of sugars into and from the phloem (Doman and Geiger, 1979; Cakmak et al., 1994; Dreyer et al., 2017). Furthermore, various processes have been described that could explain a sink limitation of tree growth under K deficiency: it could directly impact organ growth by stoichiometric limitation through its role as an enzyme co-factor or its effect on cell turgor pressure (Battie-Laclau et al., 2013), which is necessary for tissue expansion (Lockhart, 1965; Muller et al., 2011; Pantin et al., 2012). The present paper focuses on the simulation of the relationship between soil K availability and the NPP of all the main tree organs.

The objectives of the present study were to understand the following:

1. the impact of K deficiency on the NPP of a Brazilian eucalypt stand representative of large areas of commercial eucalypt plantations

2. whether the influence of K availability on GPP is sufficient to explain the differences in wood productivity between K-fertilised and K-deficient stands

3. the link between C partitioning and K availability in the plant–soil system

4. which parts of the K cycle are the most critical to simulate accurately the consequences of K limitation on wood productivity

5. if K:C stoichiometry can contribute to explaining the observed patterns of organ NPP.

We have built on the results related to the K–GPP relationship previously obtained. To this end, the C and K allocation schemes of CASTANEA-MAESPA-K were adapted to Brazilian eucalypt stands. The model was evaluated against data obtained in a K fertilisation–omission experiment (Cornut et al., 2023). To test the hypothesis of sink limitation of wood NPP, the parsimony principle was applied. In other words, we mean that we did not introduce more processes in the model if those that were already included could explain the observed patterns. Since the model did not include any explicit sink limitation process, if it was able to replicate observed productivity patterns, that meant that C-source limitation was sufficient to explain wood productivity limitation by K deficiency.

2 Methods

2.1 Study site

A split-plot fertilisation trial was installed at the Itatinga experimental station (23^∘02^′49^′′ S and 48^∘38^′17^′′ W; 860 m a.s.l; University of São Paulo-ESALQ). The precipitation was, on average, 1430 mm yr⁻¹, with a drier season between June and September, and the mean annual temperature was 19.3°C. The trial was established on June 2010 for 6 years. The planted clone was a fast growing Eucalyptus grandis. The experimental design was described in detail in Battie‐Laclau et al. (2014). Six treatments (three fertilisation regimes crossed with two water regimes) were applied in three blocks. In the present study, we focus on the +K and oK treatments with an undisturbed rainfall regime: a non-limiting K fertilisation (+K) with 17.55 gK m⁻² applied as KCl; 3.3 gP m⁻² and 200 g m⁻² of dolomitic lime and trace elements at planting and 12 gN m⁻² at 3 months of age; and an omission treatment (oK) where all fertilisation was applied as in the +K treatment, except the KCl.

The concentrations of N, P, and K (as well as of Ca and Mg) in the organs (leaves, trunks, branches, and roots) were measured at an annual time step in eight individual trees of each fertilisation treatment and were upscaled to the whole stand using allometric relationships. Biomass and nutrient contents were calculated (using upscaling) from inventories, biomass, and nutrient concentration measurements performed in each fertilisation treatment at years 1, 2, 3, 4, 5, and 6 after planting. Atmospheric deposition (0.55 gK m⁻² yr⁻¹) was measured in a nearby experiment from Laclau et al. (2010).

2.2 CASTANEA-MAESPA-K model

CASTANEA-MAESPA-K (Cornut et al., 2023) is a coupled C–H₂O–K mechanistic model that simulates forest growth at the stand level. The original CASTANEA model (Dufrêne et al., 2005), generally used to simulate temperate forest stands, was adapted to tropical eucalypt plantations. It was merged with the MAESPA model (Duursma and Medlyn, 2012; Christina et al., 2017). since it does not simulate the water–plant–atmosphere hydraulic continuum natively and since the assumption of a fixed root depth made by CASTANEA did not hold true in the studied system (Christina et al., 2011). More information on the construction of the CASTANEA-MAESPA-K model and a overview schematic can be found in Part 1 (Cornut et al., 2023) of this two-part paper.

Here, we focused on NPP (for GPP, see Cornut et al., 2022), which led to us concentrating our efforts on the C and K allocation models (Fig. 1). All the sub-models (leaf cohorts, external K cycling, uptake) described in Part 1 (Cornut et al., 2022) were re-used, albeit sometimes with different parameter values (Table S1–S2 in the Supplement), since the simulated experimental sites were different. While used implicitly in the modelling work of the companion paper (to be able to simulate K fluxes), the new modules simulating the carbon allocation, the availability of K for organ growth, and the remobilisation of K from organs developed for this work are presented below.

Figure 1 Schematic representation of K balance in the tree and its link with the allocation model. Purple boxes are K state variables, while purple arrows are K fluxes. Dashed purple arrows are remobilisation fluxes. The K uptake flux simulated with a simple Ohm's law form is represented with resistance symbols. Grey boxes represent C biomasses. Dark-grey arrows represent the allocation of NPP to the different organs. Thin dotted grey arrows represent the influence of the organ biomass on respiration.

2.3 Carbon allocation

In the CASTANEA model, carbon of assimilates was allocated to organ growth and soluble sugars (SSs) after a part of it was released to the atmosphere through maintenance respiration and growth respiration processes. SSs are not localised in the model which, for the carbon part, has no topology. So SSs are hypothesised to be part of the phloem and other tissues, indistinctly. This is different from the K content of the plant that we localised in either the organs or the semi-explicit phloem and xylem saps.

In the model, leaves had priority over other organs with regards to the allocation of carbon. This means that carbon was firstly allocated to leaves, and what was left could be allocated to the other organs. Leaf production in the model was driven by the increase in tree height (see Cornut et al., 2022, and Eq. 8 below). The parameters used here were fitted using experimental data from the +K stand. Leaf growth could, however, be limited if the demand for growth was higher than the available C in the SS compartment. The growth of all other organs was a fraction of the daily NPP minus leaf biomass production (Fig. 1). The allocation coefficients of each organ, except those of leaves (i.e. fine roots, coarse roots, woody organs), were calculated at a daily time step and were the result of NPP and allometric relationships among organs.

2.3.1 Soluble sugars

The optimal biomass of whole-tree SSs ($B_{\mathrm{SS}}^{\max }$) was a function of woody organ biomass (Eq. S1 in the Supplement). The allocation coefficient to SSs with a target value corresponding to this optimal biomass was calculated as follows:

$$\begin{array}{}\text{(1)}& G_{\mathrm{SS}}=max(\mathrm{0},G_{\mathrm{SS}}^{\max }\times min(\mathrm{1},\left(\mathrm{1}-{\displaystyle \frac{B_{\mathrm{SS}}}{B_{\mathrm{SS}}^{\max }}}\right)\times {\displaystyle \frac{\mathrm{1}}{p_{\mathrm{SS}}}}\left)\right),\end{array}$$

where G_SS (unitless) is the allocation coefficient of SSs, B_SS (gC m⁻²) is the biomass of SSs, $G_{\mathrm{SS}}^{\max }$ (unitless) is the maximum allocation coefficient of SSs, $B_{\mathrm{SS}}^{\max }$ (gC m⁻²) is the optimal biomass of SSs, and p_SS is the sensitivity of the response of the allocation coefficient to soluble-sugar deficiency. The values of the parameters of this function were chosen to ensure realistic values of SSs. Measurement-based parameterisation was impossible since no whole-tree SS measurements were conducted in eucalypts at these sites.

2.3.2 Roots

In accordance with experiment evidence in the fertilisation experiments, the allocation coefficient of coarse roots, G_CR (unitless), was set to a constant,

$$\begin{array}{}\text{(2)}& G_{\mathrm{CR}}=\mathrm{0.01}.\end{array}$$

The allocation coefficient of carbon to fine roots, G_FR (unitless), was similar to that of Marsden et al. (2013). It used a target function, where the fine-root target biomass was a function of leaf area index:

$$\begin{array}{}\text{(3)}& \begin{array}{rl}{G}_{\mathrm{FR}}=& min(\mathrm{1}-G_{\mathrm{CR}}-G_{\mathrm{SS}},max(\mathrm{0},\mathrm{0.5}\times \\ & min(\mathrm{1},\left(\mathrm{1}-{\displaystyle \frac{B_{\mathrm{FR}}}{\mathit{\lambda }\times L}}\right)\times {\displaystyle \frac{\mathrm{1}}{p_{\mathrm{FR}}}})\left)\right),\end{array}\end{array}$$

where G_CR and G_SS (unitless) are the allocation coefficients of coarse roots (Eq. 2) and SSs (Eq. 1) respectively, B_FR (gC m⁻²) is the fine-root biomass, λ (gC m${}_{\mathrm{leaf}}^{-\mathrm{2}}$) is the conversion coefficient between leaf area and target root biomass, L (m${}_{\mathrm{leaf}}^{\mathrm{2}}$ m${}_{\mathrm{soil}}^{-\mathrm{2}}$) is the leaf area index of the stand, and p_FR is a sensitivity parameter.

2.3.3 Woody organs

The remaining C available after allocation of NPP to leaves, roots, and SSs was allocated to the woody organs (trunk, branches and bark):

$$\begin{array}{}\text{(4)}& G_{\mathrm{W}}=max\left(\mathrm{0},\mathrm{1}-G_{\mathrm{FR}}-G_{\mathrm{CR}}-G_{\mathrm{SS}}\right),\end{array}$$

where G_W (unitless) is the allocation coefficient to woody organs; and G_FR, G_CR, and G_SS are the allocation to fine roots (Eq. 3), coarse roots (Eq. 2), and soluble sugars (Eq. 1) respectively.

The C allocated to woody organs was then distributed to the different woody organs (trunk, branches, and bark) using allometric relationships.

Parameters of the function linking the allocation of woody NPP to branches in the function of total woody biomass were fitted on experimental data. The following equation was used in the model:

$$\begin{array}{}\text{(5)}& G_{\mathrm{Br}}=\mathrm{0.437}\times e^{-\mathrm{0.00240}\times \left(B_{\mathrm{trunk}}+B_{\mathrm{Br}}+B_{\mathrm{bark}}\right)}+\mathrm{0.102},\end{array}$$

where G_Br (unitless) is the ratio of woody NPP allocated to branch production, and B_trunk, B_Br, and B_bark are the biomasses (gC m⁻²) of the trunk, branches, and bark respectively.

Following experimental data, the allocation of woody NPP to bark was considered to be constant.

$$\begin{array}{}\text{(6)}& G_{\mathrm{Bark}}=\mathrm{0.10}\end{array}$$

All the C remaining after allocation of woody NPP to the bark and the branches was used for trunk wood production.

$$\begin{array}{}\text{(7)}& G_{\mathrm{trunk}}=\mathrm{1}-G_{\mathrm{Bark}}-G_{\mathrm{Br}}\end{array}$$

The increase in height was a function of trunk biomass. The relationship was fitted in the +K stand using biomass and inventory data (Fig. S1 in the Supplement):

$$\begin{array}{}\text{(8)}& H=\mathrm{22.67}\times {\left(\mathrm{1}.-e^{-\mathrm{2.9}\times {\mathrm{10}}^{-\mathrm{4}}\times B_{\mathrm{trunk}}}\right)}^{\mathrm{0.4989}}+\mathrm{1},\end{array}$$

where H (m) is the height of the stand, and B_trunk (gC m⁻²) is the biomass of the trunk.

2.4 Autotrophic respiration fluxes

2.4.1 Growth respiration

Organ growth was calculated at a daily time step. In the base CASTANEA model, growth respiration was calculated using the daily growth of the organ and the associated construction cost (Dufrêne et al., 2005). In this study, the construction cost of trunk wood was based on the age of the plantation following Ryan et al. (2009):

$$\begin{array}{}\text{(9)}& {\mathrm{GRC}}_{\mathrm{trunk}}=max\left(\mathrm{0.076},\mathrm{0.368}-\mathrm{0.0343}\times \mathrm{age}\right),\end{array}$$

where GRC_trunk (gC${}_{{\mathrm{CO}}_{\mathrm{2}}}$ gC⁻¹) is the growth respiration cost associated with 1 gC of trunk wood increment, and age (days) refers to the age of the plantation. Due to the lack of available data, the growth respiration for other organs used values identical to those used in Dufrêne et al. (2005).

The daily growth respiration associated with each respective organ was the following:

$$\begin{array}{}\text{(10)}& {\mathrm{GR}}_{\mathrm{organ}}={\mathrm{GRC}}_{\mathrm{organ}}\times \mathrm{\Delta }{B}_{\mathrm{organ}}\end{array}$$

, where GR_organ (gC${}_{{\mathrm{CO}}_{\mathrm{2}}}$ m⁻² d⁻¹) is the growth respiration of the organ, GRC_organ (gC${}_{{\mathrm{CO}}_{\mathrm{2}}}$ gC⁻¹) is the growth respiration cost of the organ, and ΔB_organ (gC m⁻² d⁻¹) is the daily growth in the biomass of the organ.

2.4.2 Maintenance respiration

The hourly maintenance respiration for all organs (except leaves, which had a respiration rate based on their vertical position in the canopy; Christina et al., 2015) was a function of their respective respiration rate per nitrogen unit (see below), nitrogen content, and surface temperature (Dufrêne et al., 2005).

Firstly the respiration rates per unit biomass were calculated. For the trunk, we used data from Ryan et al. (2009) to calculate the maintenance respiration rate per unit biomass as a function of trunk biomass:

$$\begin{array}{}\text{(11)}& {\mathrm{MRN}}_{\mathrm{trunk}}=max\left(R_{\mathrm{1}}^{\mathrm{trunk}},R_{\mathrm{2}}^{\mathrm{trunk}}-R_{\mathrm{3}}^{\mathrm{trunk}}\times B_{\mathrm{trunk}}\right),\end{array}$$

where MRN_trunk (gC gN⁻¹ h⁻¹) is the respiration rate of the trunk per nitrogen unit at a reference temperature of 25 ^∘C, B_trunk is the trunk biomass in gC m⁻², $R_{\mathrm{1}}^{\mathrm{trunk}}$ (gC gN⁻¹ h⁻¹) is the minimum rate of respiration per unit nitrogen, $R_{\mathrm{2}}^{\mathrm{trunk}}$ (gC.gN⁻¹ h⁻¹) is the maximum rate of respiration per unit nitrogen, and $R_{\mathrm{3}}^{\mathrm{trunk}}$ (gN⁻¹ h⁻¹) is the slope between trunk biomass and respiration rate per nitrogen unit. Values for these parameters are reported in Table S2 in the Supplement. The respiration rate per nitrogen unit of the roots (MRN_R) and branches (MRN_Br) was assumed to be equal to twice the MRN_trunk.

Assuming a C content of 50 %, the N content of roots was fixed at 0.0050 gN gC⁻¹, since no variation was visible in the experimental data. The N content of branches and trunks was a decreasing exponential function of their respective biomasses that was calibrated using experimental data as follows (Fig. S2 in the Supplement):

$$\begin{array}{}\text{(12)}& \begin{array}{rl}{N}_{\mathrm{trunk}}& =min(\mathrm{9.50}\times {\mathrm{10}}^{-\mathrm{3}},\mathrm{9.92}\times {\mathrm{10}}^{-\mathrm{3}}\\ & \times e^{-\mathrm{4.13}\times {\mathrm{10}}^{-\mathrm{3}}\times B_{\mathrm{trunk}}}+\mathrm{1.61}\times {\mathrm{10}}^{-\mathrm{3}})\\ N_{\mathrm{branches}}& =min(\mathrm{1.10}\times {\mathrm{10}}^{-\mathrm{2}},\mathrm{5.12}\times {\mathrm{10}}^{-\mathrm{3}}\\ & \times e^{-\mathrm{9.66}\times {\mathrm{10}}^{-\mathrm{3}}\times B_{\mathrm{branches}}}+\mathrm{5.21}\times {\mathrm{10}}^{-\mathrm{3}}),\end{array}\end{array}$$

where N_trunk and N_branches are the respective N concentrations (gN gC⁻¹) of the trunk and the branches, and B_trunk and B_branches are the respective biomasses (gC m⁻²) of the branches and the trunk. The N content of leaves and bark was not simulated, since N content did not influence the respiration of leaves in our model and since bark had no maintenance respiration.

Following Dufrêne et al. (2005), the maintenance respiration of an organ (except leaves, see below) was

$$\begin{array}{}\text{(13)}& {\mathrm{RM}}_{\mathrm{organ}}=B_{\mathrm{organ}}\times {\mathrm{MRN}}_{\mathrm{organ}}\times N_{\mathrm{organ}}\times Q_{\mathrm{10}}^{\left(T_{\mathrm{org}}-T_{\mathrm{MR}}\right)/\mathrm{10}},\end{array}$$

where RM_organ (gC h⁻¹ m⁻²) is the respiration rate of the organ, B_organ (gC m⁻²) is the biomass of the organ, MRN_organ is the respiration rate per unit nitrogen (gC gN⁻¹ h⁻¹, Eq. 11), N_organ (gN gC⁻¹) is the N concentration of the organ, and Q₁₀ is the exponential relationship between the respiration rates and temperature.

The maintenance respiration of leaves was their dark respiration (inhibited during the day), R_d, and the values measured in K-fertilised trees in a nearby site (Eucflux) were used (Fig. S1 in Christina et al., 2017).

2.5 Organ turnover

With the exception of the coarse roots and the trunk, the organs (branches, bark, fine roots) were subject to turnover. Branches, bark, fine roots, and leaves each had lifespans.

In simulations, the theoretical (independent of K limitation effects) leaf lifespan (LLS) was considered to be constant throughout the rotation. This option was chosen because we were unable to mechanistically model the observed variations in leaf lifespan. The realised lifespan of leaves was influenced by their K status (leaves fell when their K concentration was below a certain threshold; see Sect. 2.5.3 in Cornut et al., 2022). We assumed that the lifespans of bark (BarkLS, Fig. S3b in the Supplement), branches (BrLS, Fig. S3a in the Supplement), and fine roots (FRLS, Lambais et al. (2017)) depended on neither tree age nor the nutritional status of the organs.

The necromass of most organs was added to the litter pool. On the other hand, dead branches were added to the dead-branch pool. The dead-branch pool represented the branches that stay attached to the tree after senescence. The dead-branch pool had a specific turnover rate.

Resorption of K took place during the senescence of leaves (Cornut et al., 2022) and branches (Fig. S4a in the Supplement). K remobilised from these two organs was added to the phloem sap K pool. While the resorption rate for leaves was dependent on the nutritional status of the tree and their theoretical lifespan, it was fixed for branches.

2.6 K allocation

The allocation of K to organs was a function of the optimal K concentration of newly formed organ tissue, organ NPP, and K availability in the tree. Firstly, organ NPP was calculated by allocating part of GPP to the organs after subtracting the respiration. Then, the total K quantity required was calculated by multiplying the growth of each organ by its optimal K concentration (the concentration of newly formed organs in the fully fertilised stand). If the quantity of available K in the phloem sap was inferior to the demand, K allocation to the organs was limited without affecting C allocation to the organ's growth. This was equivalent to flexible stoichiometry in other models.

$$\begin{array}{}\text{(14)}& \begin{array}{rl}{K}_{\mathrm{NPP}}& ={\mathrm{NPP}}_{\mathrm{org}}\times \left(\right[K{]}_{\mathrm{trunk}}^{\mathrm{opti}}\times G_{\mathrm{trunk}}+[K{]}_{\mathrm{Br}}^{\mathrm{opti}}\times G_{\mathrm{Br}}\\ & +[K{]}_{\mathrm{Bark}}^{\mathrm{opti}}\times G_{\mathrm{Bark}}+[K{]}_{\mathrm{CR}}^{\mathrm{opti}}\times G_{\mathrm{CR}}\\ & +[K{]}_{\mathrm{FR}}^{\mathrm{opti}}\times G_{\mathrm{FR}}),\end{array}\end{array}$$

where K_NPP (gK m⁻² d⁻¹) is the quantity of K necessary for optimal stoichiometry of newly formed organ biomass, NPP_org (gC m⁻² d⁻¹) is the daily net primary production minus the allocation to the leaves, K${}_{\mathrm{org}}^{\mathrm{opti}}$ (gK gC⁻¹) is the optimal concentration of the considered organ (Eqs. 17–19), and G_org is the allocation coefficient of that organ (Eqs. 2–4).

The quantity of available K was a function of K content in the phloem sap and the minimal quantity of K in the phloem sap (Cornut et al. 2022):

$$\begin{array}{}\text{(15)}& K_{\mathrm{available}}={\displaystyle \frac{K_{\mathrm{phloem}}-K_{\mathrm{phloem}}^{\min }}{\mathrm{\Delta }t}},\end{array}$$

where K_available (gK m⁻² d⁻¹) is the amount of K available for organ growth, K_phloem (gK m⁻²) is the amount of K in the phloem sap, K${}_{\mathrm{phloem}}^{\min }$ (gK m⁻²) is the minimal amount of K in the phloem sap (Cornut et al., 2022), and Δt (days) is the time step (1 d in our simulations).

The limitation of K allocation to newly formed organ nutrient content was simply the ratio between available K and K demand:

$$\begin{array}{}\text{(16)}& L_{\mathrm{K}}={\displaystyle \frac{min(K_{\mathrm{available}},K_{\mathrm{NPP}})}{K_{\mathrm{NPP}}+K_{\mathrm{Leaf}}^{\mathrm{Demand}}}},\end{array}$$

where K_available is from Eq. (15), K_NPP (gK m⁻² d⁻¹) is the amount of K needed for optimal stoichiometry of newly formed woody organs (Eq. 14), and $K_{\mathrm{Leaf}}^{\mathrm{Demand}}$ (gK m⁻² d⁻¹) is the leaf growth K demand (Cornut et al., 2023). The cycle of K in the leaves is described in Cornut et al. (2022), since it is an integral part of the canopy cohort model.

2.6.1 Trunk wood

Due to the continuous phenology of tropical eucalypt trees, K dynamics in trunk wood were simulated through daily cohorts of trunk wood. This was also the preferred modelling choice, since it provided a mechanistic explanation for trunk wood remobilisation throughout trunk wood ageing (that we considered to be leaching of K from sapwood into the xylem sap). The trunk NPP was allocated daily to a newly formed cohort of trunk wood. In parallel, the optimal K concentration of newly formed trunk wood ([K]${}_{\mathrm{trunk}}^{\mathrm{opti}}$) was constant and equal to the maximum trunk wood concentration measured by the fertilisation experiment (Fig. S5 in the Supplement). The realised concentration of newly formed trunk wood was a function of the optimal K concentration of newly formed trunk wood and the strength of K supply limitation; i.e. the K content of cohort i at its creation was the following:

$$\begin{array}{}\text{(17)}& K_{\mathrm{phloem}\to \mathrm{trunk}}^i=[K{]}_{\mathrm{trunk}}^{\mathrm{opti}}\times L_{\mathrm{K}}\times {\mathrm{NPP}}_{\mathrm{trunk}},\end{array}$$

where $K_{\mathrm{phloem}\to \mathrm{trunk}}^i$ (gK m⁻² d⁻¹) is the flux of K to the cohort i at the time of its creation, [K]${}_{\mathrm{trunk}}^{\mathrm{opti}}$ (gK gC⁻¹) is the optimal concentration of newly formed trunk wood, Lim${}_{\mathrm{org}}^{\mathrm{K}}$ is the limitation by K availability (Eq. 16), and NPP_trunk (gC m⁻² d⁻¹) is the daily trunk wood increment.

2.6.2 Branches

The optimal K concentration of newly formed branches was a function of the branch biomass,

$$\begin{array}{}\text{(18)}& [K{]}_{\mathrm{Br}}^{\mathrm{opti}}=\mathrm{1.41}\times {\mathrm{10}}^{-\mathrm{2}}\times e^{-\mathrm{8.54}\times {\mathrm{10}}^{-\mathrm{3}}\times B_{\mathrm{Br}}}+\mathrm{1.52}\times {\mathrm{10}}^{-\mathrm{3}},\end{array}$$

where [K]${}_{\mathrm{Br}}^{\mathrm{opti}}$ is the K concentration of the newly formed branches in gK gC⁻¹, and B_Br is the branch biomass in gC m⁻². This decreasing function was fitted on experimental nutrient content and biomass data collected in the fertilised plots.

2.6.3 Bark

The optimal K concentration of newly formed bark was a function of the bark biomass,

$$\begin{array}{}\text{(19)}& [K{]}_{\mathrm{Bark}}^{\mathrm{opti}}=\mathrm{4.08}\times {\mathrm{10}}^{-\mathrm{3}}\times e^{-\mathrm{5.68}\times {\mathrm{10}}^{-\mathrm{3}}\times B_{\mathrm{Bark}}}+\mathrm{2.95}\times {\mathrm{10}}^{-\mathrm{3}},\end{array}$$

where [K]${}_{\mathrm{Bark}}^{\mathrm{opti}}$ is the K concentration of the newly formed bark in gK gC⁻¹, and B_Bark is the bark biomass in gC m⁻². The parameters for this function were fitted on experimental data collected in the fertilised plots. No remobilisation was considered for bark, since there no measurements were available.

2.6.4 Roots

The optimal K concentration (in gK gC⁻¹) of coarse ([K]${}_{\mathrm{CR}}^{\mathrm{opti}}$) and fine roots ([K]${}_{\mathrm{FR}}^{\mathrm{opti}}$) was a fixed value independent of tree age or biomass. Due to the absence of data regarding this process, the model did not simulate remobilisation from roots. The K content of dead fine roots was added to the K litter pool, which in turn leached into the soil available K and could be taken up by other living roots.

2.7 K remobilisation and turnover

2.7.1 Branches

Before falling as litter to the ground, branches that die stay attached some months to the tree. Whenever dead branches are mentioned, we mean dead branches still attached to the trunk. During branch death, part of the branch K was remobilised into K_xylem. The remobilisation efficiency (R_KBranch, Eq. 20) was considered to be constant and was calculated using the difference between the dead-branch and live-branch K concentration in experimental data. The value of remobilisation rate used in the model was calculated based on the K remobilisation rate of branches during the first year of the rotation to avoid including the impact of leaching from dead branches (Fig. S4 in the Supplement). Indeed, we postulated that the dead branches present in the first year of the rotation were younger than those present at later stages, reducing the potential for leaching flux and the impact on branch K remobilisation calculations. Leaching fluxes from dead branches remaining attached to the tree trunk were assumed to be linearly related to the duration of time a dead branch remained on the tree.

2.7.2 Trunk

Since the K concentration of the total trunk tissue decreases with trunk biomass, it was necessary to implement trunk wood K remobilisation in the model. A model where the remobilisation rate was dependent on trunk wood production was the best suited for this task, since it showed the best fit with experimental data when compared to a model where remobilisation was independent of trunk wood production. The following equation was used:

$$\begin{array}{}\text{(20)}& K_{\mathrm{trunk}\to \mathrm{xylem}}^i=K_{\mathrm{trunk}}^i\times T_{\mathrm{Ktrunk}}\times {\mathrm{NPP}}_{\mathrm{trunk}},\end{array}$$

where K${}_{\mathrm{trunk}\to \mathrm{xylem}}^i$ (gK m⁻² d⁻¹) is the remobilisation of K from the cohort i, K${}_{\mathrm{trunk}}^i$ (gK m⁻²) is the K mineralomass of the trunk wood cohort i, T_Ktrunk (gC⁻¹ m²) is the remobilisation rate per unit of trunk wood production, and NPP_trunk (gC m⁻² d⁻¹) is the daily trunk increment. Remobilised K was allocated to K_xylem. If the K concentration of the cohort ([K]${}_{\mathrm{trunk}}^i$) was lower than a threshold value [K]${}_{\mathrm{trunk}}^{\min }$ (gK gC⁻¹), there was no remobilisation (K${}_{\mathrm{trunk}\to \mathrm{xylem}}^i=\mathrm{0}$). The threshold, [K]${}_{\mathrm{trunk}}^{\min }$, was determined from the minimum asymptote of the relationship between trunk wood biomass and trunk wood K concentration in the fertilisation experiment (Fig. S5 in the Supplement). This measured value was assumed to be the minimum concentration of a cohort, since it was assumed that, at a high-enough wood biomass, the proportion of K associated with newly formed wood to the total wood K content was negligible (Augusto et al., 2000). This meant that the measured concentration of trunk wood as a whole was similar to the minimum concentration of K in the trunk at high-enough trunk biomass.

2.7.3 Total remobilisation in woody organs

Total remobilisation was the flux of K from the woody organs to the xylem. It was calculated as the following:

$$\begin{array}{}\text{(21)}& K_{\mathrm{remob}}=K_{\mathrm{trunk}\to \mathrm{xylem}}+R_{\mathrm{Kbranches}}\times K_{\mathrm{branches}}^{\mathrm{mortality}}\end{array}$$

where K_remob (gK m⁻² d⁻¹) was the total remobilisation flux, K_{trunk→xylem} (gK m⁻² d⁻¹) the remobilisation rate of wood, R_Kbranches (unitless) the remobilisation rate of dying branches and K${}_{\mathrm{branches}}^{\mathrm{mortality}}$ (gK m⁻² d⁻¹) the flux of K from living branches to dead branches.

2.8 Simulations

The simulation initialisations were conducted to resemble as closely as possible the omission experiment. Simulations in the fully fertilised treatment (+K) were initialised with the same fertilisation values as the fertilised control in the experiment (i.e. 17.5 gK m⁻²), corresponding to a one-time application of fertiliser at planting. Simulations in the K omission treatment (oK) shared the same initialisation, except that the fertiliser pool was initialised with 0 gK.

To investigate the effects of a fertilisation gradient, 10 initialisation values of K pools spanning from no input (in oK) to 17.5 gK (+K) were chosen.

To test whether the fertilisation regime could have an impact on tree productivity, the following two fertilisation regimes were simulated: one where the K dose was brought all at once (as in the experiment) and one where the same K fertiliser dose was broken up into four sub-doses that were temporally spaced (equivalent to the Eucflux experiment; Cornut et al., 2022).

2.9 Analysis

To test the accuracy of the model prediction, the root-mean-square errors of the simulations' output variables were calculated using measurements at the experimental site. The mean of the three experimental blocks (there were three blocks per experimental treatment) was used. To normalise this metric and to have a relative root-mean-square error, the RMSE was divided by the measured mean throughout the rotation of the considered output variable.

To describe the response of resource use efficiency (RUE) to different levels of K availability, we used the following metrics.

2.9.1 Carbon use efficiency

The carbon use efficiencies (CUE) were calculated as the considered simulated C flux summed over the whole rotation (CUE_NPP for NPP and CUE_trunk for trunk NPP) divided by the simulated GPP flux summed over the whole rotation (De Lucia et al., 2007). It is in fact a measure of the proportion of assimilated carbon that was used for forming tissue, i.e. not re-emitted through autotrophic respiration.

2.9.2 Water use efficiency

The water use efficiencies of NPP (WUE_NPP), trunk NPP (WUE_trunk), and GPP (WUE_GPP) were calculated by dividing the total NPP, trunk NPP, and GPP respectively by the amount of transpired water during the period over which NPP, trunk NPP, and GPP were calculated (here the whole rotation).

2.9.3 Potassium use efficiency

Potassium use efficiencies of GPP (KUE_GPP), total NPP (KUE_NPP), and trunk NPP (KUE_trunk) were calculated by dividing the respective C flux by the maximum amount of K that was immobilised in the plant during the rotation. For example, in the case of (KUE_NPP),

$$\begin{array}{}\text{(22)}& {\mathrm{KUE}}_{\mathrm{NPP}}={\displaystyle \frac{{\sum }_{t=\mathrm{0}}^k{\mathrm{NPP}}_t}{K_{\mathrm{plant}}^{\max }}},\end{array}$$

where KUE_NPP (gC gK⁻¹) is the K use efficiency of total NPP, k is the number of days in the rotation (days), NPP_t (gC m⁻²) is the daily NPP of the rotation, and K${}_{\mathrm{plant}}^{\max }$ (gK m⁻²) is the maximum of K that was immobilised in the plant during the rotation (the maximum of total simulated plant K during the rotation). For calculating KUE_GPP and KUE_trunk, the numerator of the above fraction can be replaced by GPP or NPP_trunk respectively. There are many alternative ways to calculate nutrient use efficiencies in forests (Turner and Lambert, 2014). Here, we decided to use total K immobilisation instead of uptake, since circulation of K in the system was high, and we think that the maximum amount of K accumulated in standing biomass is a more relevant representation of total system K demand. Indeed, the maximum K accumulated in standing biomass is a proxy of the amount of K necessary in the system and along the rotation for the plant considering its biomass. In systems where there is less restitution to the soil, it would be equivalent to the soil nutrient uptake.

2.9.4 Fertiliser use efficiency

Fertiliser use efficiencies were computed as the difference of cumulated NPP between the simulated K omission stand (oK) and stands simulated with different K fertilisation levels, divided by the amount of K fertiliser added. This allowed us to compute the growth gain in carbon per unit of K fertiliser used as follows:

$$\begin{array}{}\text{(23)}& {\mathrm{FUE}}_{\mathrm{NPP}}^f={\displaystyle \frac{{\sum }_{i=\mathrm{0}}^k\left({\mathrm{NPP}}_i^f-{\mathrm{NPP}}_i^{\mathrm{oK}}\right)}{{\sum }_{i=\mathrm{0}}^k{K}_{\mathrm{fertiliser},i}^{\mathrm{added}}}},\end{array}$$

where FUE${}_{\mathrm{NPP}}^f$ (gC gK⁻¹) is the fertiliser use efficiency of NPP for a given level of fertilisation, k is the number of days in the rotation (days), NPP^f (gC m⁻² d⁻¹) is the daily NPP of the currently considered stand, NPP^oK (gC m⁻² d⁻¹) is the NPP of the K omission stand, and K${}_{\mathrm{fertiliser},i}^{\mathrm{added}}$ (gK m⁻² d⁻¹) is the amount of K fertiliser that was added at day i in the considered stand. To obtain FUE_GPP or FUE_trunk, this relationship can be applied to either GPP or NPP_trunk respectively.

3 Results

3.1 Prediction of changes in NPP caused by GPP

The model was capable of replicating most of the NPP and biomass differences between the +K and oK stands (Fig. 2). In the +K stand, the 5-year yearly averaged GPP, NPP, and trunk NPP (NPP_trunk) fluxes were respectively 3966, 1990, and 1159 gC m⁻² yr⁻¹. In the oK stand, they were respectively 1781, 715, and 414 gC m⁻² yr⁻¹. The GPP, NPP, and NPP_trunk were respectively 55 %, 64 %, and 70 % lower in the oK stand compared to in the +K stand. The reduction in GPP was in line with what was simulated at the Eucflux site (Table 3 of Cornut et al., 2022) and GPP estimations using the total below-ground carbon allocation (TBCA) method applied to C stocks and C fluxes measured throughout the rotation in our experiment (not shown; Giardina and Ryan, 2002). The reduction of NPP_trunk was comparable to data (Fig. 2b). The same could be said for the bark (Fig. 2d) and the branches (Fig. 2f). This led to simulated biomasses in line with measurements in the +K and oK stands for branches (Fig. 2e), bark (Fig. 2c), and the trunk (Fig. 2a). In the +K stand, the RMSEs (and normalised RMSEs in parentheses) of the simulated trunk, branches, bark, leaves, and total aboveground biomass were respectively 385 gC m⁻² (12 %), 50 gC m⁻² (17 %), 34 gC m⁻² (9 %), 48 gC m⁻² (19 %), and 363 gC m⁻² (9 %). In the oK stand, they were respectively 155 gC m⁻² (16 %), 48 gC m⁻² (37 %), 25 gC m⁻² (17 %), 31 gC m⁻² (26 %), and 193 gC m⁻² (15 %). The errors of the simulations relative to measured values of total aboveground biomass at month 59 after planting were an underestimation of 24 gC m⁻² (−1 %) in the +K stand and an overestimation of 161 gC m⁻² (6 %) in the oK stand. This overestimation of aboveground biomass was concurrent to an underestimation of root biomass (Fig. 2i).

Figure 2 Measured and simulated biomass of the different organs in two contrasted K availability scenarios (no K fertilisation, oK, and full K fertilisation, +K) in (a) trunk, (c) bark, (e) branches, (g) leaves, and (i) fine roots. In the right column are the respective measured and simulated annual NPP of each of the following organs: (b) trunk, (d) bark, (f) branches, (h) leaves, and (j) fine roots.

3.2 Consequences of K addition for C allocation patterns within trees

The model allowed the study of allocation in the trees under different K fertilisation regimes. The simulated allocation patterns did not differ greatly between the fertilised and omission stands (Fig. 3a). However, simulated CUE_NPP was reduced by 23 % in the omission stand (Fig. S6b in the Supplement; 0.40 vs. 0.52). The ratio of wood productivity to GPP (CUE_trunk) was reduced by the same proportion (21 %), showing that the reduction of NPP_trunk followed the same dynamic as the total NPP. Moreover, the difference in CUE_NPP between the two fertilisation treatments increased throughout the rotation (a 2 % difference the first year and a 17 % difference the fifth year). The trend was similar for CUE_trunk. The difference in CUE_NPP between the two treatments was mainly due to a relative increase of maintenance respiration in the oK stand, where it accounted for 48 % of the GPP compared to 34 % in the +K stand. For CUE_trunk, this was further amplified by leaf NPP representing 13 % of GPP in oK compared to 7 % in +K.

The model was also capable of simulating the response of the different carbon fluxes, along the stand rotation, to a gradient of initial K fertilisation (Fig. 3b). It showed that the responses of GPP, NPP, and wood productivity to fertilisation all saturated at around 11 gK m⁻² for a 5-year rotation. The simulated carbon fluxes did not show any sensitivity to the fertilisation application regime (one- or four-time application; Fig. S6a in the Supplement). The simulations conducted with the one-time application at planting compared with the same amount of K split into four applications at months 0, 3, 10, and 20 after planting showed little to no difference in GPP, NPP, and NPP_trunk (Fig. S6a in the Supplement). The wood productivity was similar in all fertilisation treatments in the first year of the rotation (Fig. 3). This result was in contrast to experimental data that show that the relative difference in organ NPP appears early in the rotation (Laclau et al., 2009). While the response of CUE_NPP resembled a linear function before it saturated at a fertilisation of 11 gK m⁻², the response of CUE_trunk followed a non-linear response (Fig. S6b in the Supplement) by increasing from 0 to 2 gK m⁻² of fertilisation added as KCl fertiliser, saturating between 2 and 4 gK m⁻² and increasing linearly between 4 and 11 gK m⁻².

Figure 3 (a) The simulated allocation of assimilated carbon to the different organs and respiration fluxes in two contrasting fertilisation conditions: +K and oK (hatched). (b) The response of simulated carbon fluxes, cumulated over 5 years, to different fertilisation levels. Total GPP is in green, total NPP is in black, and trunk NPP is in brown.

3.3 K cycling in the trees

The model showed that the main sinks of K in both the +K and oK stands were located in the woody organs (trunk, bark, and branches). Despite the remobilisation of K in the trunk, the quantity of K immobilised in the trunk increased linearly with time in both treatments (Fig. 4), thus constituting an important K sink.

The theoretical minimum concentration of K in the xylem sap (assuming no recirculating K) of our trees was calculated by dividing the daily simulated flux of K that circulated in the xylem sap (uptake, wood, and branch remobilisation) by the simulated transpiration flux of each day. The mean simulated minimum xylem sap K concentration over the course of a rotation was 0.30 mM (0.012 gK L⁻¹) in the fully fertilised stand and 0.11 mM (0.004 gK L⁻¹) in the K omission stand. When including the K content of the phloem sap and leaf resorption (which means the total circulating K in the tree) in this calculation, the values were 1.66 mM (0.065 gK L⁻¹) and 0.46 mM (0.018 gK L⁻¹) respectively.

The simulation of internal and external K fluxes in the system (Table 1) showed that, in the fully fertilised and K omission stands, wood remobilisation represented the most important flux of K. Implementing this in the model, this process of K remobilisation from wood increased the model accuracy substantially (not shown here) by buffering the amount of K available for organ growth. When added to branch and leaf resorption, the total amount of K remobilised represented 1.8 times the K uptake in the fertilised stand versus 1.4 times the K uptake in the omission stand. In the simulated oK stand, K uptake was very similar to the sum of the litterfall, leaching, and atmospheric deposition fluxes in all but the first year of the rotation (Table 1). The deposition flux represented more than 50 % of the uptake flux in these K-deficient conditions. Moreover, in the K omission stand, an increase of the simulated weathering flux from 0 to 0.3 gK m⁻² yr⁻¹ (in the range of possible values; see Cornut et al., 2021) led to an increase of 23 %, 28 %, and 30 % of the rotation-cumulated GPP, NPP, and NPP_trunk respectively. This showed that small differences in K input can lead to big differences in outcome for wood productivity.

The difference in the K flux of canopy leaching between the two simulated stands (much lower in the omission stand, Table 1) was in line with results obtained on eucalypt plantations at K-rich and K-deficient sites (Laclau et al., 2010). This was the result of lower leaf K concentration and supports the validity of the leaching model used here (Cornut et al., 2022).

Figure 4 Simulated K mass of all the compartments containing K. Two contrasted K availability scenarios are displayed: (a) high fertilisation (+K in the main text) and (b) no K fertilisation (oK in the main text). Note the differences in the y-axis scale. The large amount of K present in the litter compartment (solid grey line) at the beginning of the rotation was measured in the plots at planting (forest floor + aboveground harvest residues).

Table 1 Simulated yearly fluxes of K (in gK m⁻² yr⁻¹) in a fully fertilised treatment (+K) and in a K omission treatment (oK). This table contains both internal fluxes (wood, branch, and leaf remobilisations) and exchange fluxes (uptake, litterfall, and canopy leaching). The constant atmospheric deposition flux (0.5 gK m⁻² yr⁻¹) is not shown.

3.4 Water and potassium use efficiencies

Omission of K fertiliser decreased stand transpiration by 51 %. The reduction (13 %) of simulated GPP water use efficiency (WUE_GPP) between the +K (0.0035 gC_GPP gH₂O⁻¹) and oK (0.0031 gC_GPP gH₂O⁻¹) stands was comparable to results obtained at the Eucflux site (Cornut et al., 2022). On the other hand, simulated WUE_NPP showed a much stronger response to K deficiency. The reduction was on the order of 33 %, with WUE_NPP at 0.0018 and 0.0012 gC gH₂O⁻¹ in the +K stand and oK stand respectively. The simulated water use efficiency of trunk wood was also reduced by 32 % in oK relative to +K. The WUE_trunk values in the +K and oK stands were respectively 0.0011 and 0.0007 gC gH₂O⁻¹.

In the following paragraph, potassium use efficiency (KUE) is understood as the ratio of accumulated carbon flux at the end of the rotation to the maximum value of K immobilised in the tree. The simulated KUE_GPP values were 1281 and 1994 gC gK${}_{\mathrm{plant}}^{-\mathrm{1}}$ in the +K and oK stands respectively. In contrast, the simulated KUE_NPP and KUE_trunk only increased by 19 % (owing to decreased CUE in the oK stand) between the +K and oK stands. The simulated KUE_trunk values were 387 and 462 gC gK${}_{\mathrm{plant}}^{-\mathrm{1}}$ in the +K and oK stands respectively (656 and 784 respectively for KUE_NPP).

FUE generally decreased with increasing fertilisation (Fig. S7 in the Supplement). However, a two-slope relationship was apparent: at low levels of fertilisation, increases in fertilisation led to strong increases in NPP. However, for high amounts of fertiliser applied, wood production per unit of K declined. At intermediate levels of fertilisation (between 6 to 10 gK m²), the FUE was almost constant or slightly increased. However, at higher levels of initial fertilisation, the FUE linearly decreased (consequence of a stable NPP with linearly increasing fertilisation, Fig. 3b).

3.5 Stoichiometry of organs

The new model of trunk wood growth and K remobilisation was validated by the simulated concentrations that were in line with experimental measurements in both the fully fertilised and the K omission treatments (Fig. 5a) without the need for additional forcing. However, the simulated concentration of K in branches (Fig. 5b) was overestimated in the oK compared to measurements. The model was unable to replicate the difference in terms of K concentration in branches between the +K and oK stands at the beginning of the rotation. Branches, while having a small biomass compared to the trunk, were an important K stock in the simulated oK stand (Fig. 4b), and it was the most important K stock at the beginning of the oK rotation. Measurements of K concentration in branches were highly variable at lower branch biomass (Fig. 5b). The mean K concentrations in the total tree biomass were 0.0048 gK gDM⁻¹ in the +K stand and 0.0030 gK gDM⁻¹ in the oK stand. This corresponded to a decrease of 36 % of the K concentration in the oK stand relative to the +K stand. This revealed that total plant stoichiometric flexibility was high in the model, in accordance with measurements.

Figure 5 Simulated concentrations of K in stem wood (a) and branches (b) in two contrasting K availability scenarios compared to measurements conducted in the fertilisation experiment at Itatinga.

4 Discussion

4.1 GPP limitation of NPP and partitioning of photosynthates

The CASTANEA-MAESPA-K model was largely successful in reproducing the limitation of wood productivity induced by K deficiency (Fig. 2). This makes it the first mechanistic model to simulate the interaction between the K cycle and forest NPP. Combined with the fact that partitioning in the model was not directly impacted by K availability, there was no mechanism through which K directly impacted carbon allocation to wood, which suggests that the limitation of wood productivity in the absence of K fertilisation was mainly due to GPP limitation. However, the consequences of K deficiency were higher for NPP than for GPP. This was due to a decrease in CUE_NPP at low levels of fertilisation. This is similar to the general trend for world forests, which is a decrease in forests' biomass-production–GPP ratio (a proxy for CUE) with decreasing fertility (Vicca et al., 2012). Here, the partitioning of GPP to the different organs was not strongly affected by K availability (except for leaves, going from 7 % of NPP in +K to 13 % of NPP in oK). The reduced CUE_NPP was mainly the result of an increased autotrophic respiration in proportion to the biomass of organs mirroring the increase in the ratio of ecosystem respiration to GPP in nutrient-poor forests (Fernández-Martínez et al., 2014). While the relationship between GPP and net ecosystem productivity was not significant for low-fertility forest sites in this meta-analysis, this was not the case for GPP and NPP in our study.

Only simulated fine-root biomasses were underestimated in both the +K and oK stands (Fig. 2i). This suggests that simulating a root target biomass that is a function of leaf area (Marsden et al., 2013) might not be appropriate in this instance (since the shape of the root biomass curve during the rotation is qualitatively different). One other cause could be a misestimation of root lifespan, since measurements of fine-root turnover have yielded a wide range of values (Jourdan et al., 2008; Lambais et al., 2017). The estimation of fine-root turnover had also been a constraint in simulating N mineralisation rates in Australian eucalypt stands with the G'DAY model (Corbeels et al., 2005). The values of GPP that were simulated here (a mean of 3986 gC m⁻² yr⁻¹ in +K and 1709 gC m⁻² yr⁻¹ in oK) were in the high range of values expected for terrestrial ecosystems (Baldocchi and Penuelas, 2019; Luyssaert et al., 2007), especially considering these values were a 6-year mean that included the first year after planting. They were similar to values estimated using TBCA measurements in highly productive eucalypt plantations (Ryan et al., 2004, 2010) and in clonal experiments at the Eucflux site (unpublished data; see Cornut et al., 2022 for the site description).

The K fertilisation level at which K limitation was totally alleviated in simulations (Fig. 3b) was similar to the 12 gK m⁻² of fertilisation that is commonly added to commercial eucalypt plantations in Brazil (Cornut et al., 2021).

4.2 Water and potassium use efficiency is affected by K availability

The simulated RUEs of NPP and NPP_trunk were strongly affected by K availability (modelled here as a change in fertilisation levels). Variations of WUE_trunk are in line with experimental results that showed a decrease in WUE_trunk in the K omission stand at the Itatinga site (Battie-Laclau et al., 2016). While simulated WUE_trunk was 33 % lower in oK than in +K, measurements showed a decrease of 37 % (Battie-Laclau et al., 2016). This confirms the relevance of a model-based approach in studying the effect of K availability on WUE_trunk. However, the model diverged from measurements in the same stands at the end of the stand rotation (from 4 to 6 years after planting) that showed that the WUE_trunk of the K omission stand was reduced by 75 % as compared to the WUE_trunk of the +K stand (Asensio et al., 2020), while CASTANEA-MAESPA-K only showed a reduction of 39 % of WUE_trunk between these two stands at this age. Since the trunk NPP simulated by the model agreed with measurements (Fig. 2b) and the model simulated a reduction in transpiration consistent with different approaches (Cornut et al., 2022), WUE_trunk in the oK stand might have been underestimated by the measurements presented in Asensio et al. (2020).

The changes in both KUE and FUE along a K fertilisation gradient showed some interesting results. While K fertilisation strongly decreased KUE for GPP, the effects on NPP and trunk NPP were weaker as a result of increases in the respective CUEs (Fig. S7 in the Supplement). FUE results demonstrate that the response of wood productivity is not a linear function of fertilisation and that, at low levels of fertilisation, small increases in K fertilisation levels produce a strong increase in NPP. Simulated NPP was maximal with the highest fertiliser amount; however, the trend is asymptotical, with more than 95 % of trunk NPP (Fig. 3b) already having reached around 10 gK m⁻², a value commonly applied in commercial eucalypt plantations managed with this soil type. Partitioning of the K fertilisation in only one or in several amounts did not changed the wood production, considering that no deep leaching occurs in these systems, in agreement with the conclusions of field studies measuring soil solution chemistry in deep soil layers (Laclau et al., 2010).

4.3 Circulation of K in the plant and stoichiometry

Our work pinpointed the importance of a plant K circulation model. The total remobilisation flux of K from the branches, the trunk, and the leaves was higher than the K uptake in the soil at all fertilisation levels. In C–N (Zaehle et al., 2010; Thum et al., 2019) or C–N–P (Goll et al., 2017) coupled models, stoichiometry of organs can be a direct limiter of organ growth. A strong effect of stoichiometry on soil organic matter decomposition has also been historically used in models of organic matter decomposition in soils (Parton et al., 1988). While such mechanisms were not considered in our modelling approach (organ K concentration is allowed to vary unconstrained), the reduction in NPP was enough to compensate for the reduction in K availability such that the simulated stoichiometries of organs do not vary more than they do in the measurements (Fig. S5 in the Supplement). The stoichiometric flexibility of trees could be higher for K than for N and P. This is consistent with observations of wood K content, which depends mainly on abiotic conditions, while P wood content depends mainly on the species (Bauters et al., 2022). Our model suggests that leaves have a higher stoichiometric flexibility than wood (not shown here), in agreement with observations in a Mediterranean forest environment (Sardans et al., 2012). However, the model failed at reproducing the patterns of stoichiometric plasticity (variability in K concentration of organs) that were observed in between the woody organs, particularly in branches (measurements show high stoichiometric flexibility in branches). This suggests that organs can differ in K homeostasis. The failure of the model to reproduce stoichiometric flexibilities in wood and bark (not shown here) could potentially have had an influence on the amount of K available for leaf expansion, thus overestimating the K limitation on canopy surface and leaf functioning. This is further exacerbated by the large amount of K stored in the bark and in the branches in the simulated K omission stand, especially when the model showed a strong K limitation in leaves (between the 10th and 20th months after planting).

The simulated amount of K immobilised in the trunk at the end of the rotation in the +K stand was 1 order of magnitude lower than what is observed in tropical forests (Bauters et al., 2022). This is due to a low total biomass in planted forests managed in short rotation compared to a in natural forest, since the simulated K concentrations in the trunk were in the range of values reported in old-growth tropical forests. The yearly increase in the amount of K stored in the trunk was in the range of observed values (see Fig. 1 in Bauters et al., 2022) for both the oK and +K stands. This suggests that tropical eucalypt plantations could be a relevant model system for certain parts of the K cycle. The quantity of K that was allocated daily to wood and remobilised gradually (see Eq. 7) acted as a buffer that prevented an overestimation of K limitation once K in the soil reached very low values. This shows that wood can act as a storage organ for K and can alleviate low uptake of K under drought (Sardans and Peñuelas, 2007; Touche et al., 2022) or in planted forests only fertilised at planting and growing in K-poor soils.

While the K trunk remobilisation fluxes were higher in our model than what was previously calculated using other methods (a mean of 7.6 gK m⁻² yr⁻¹ vs. 2.7 gK m⁻² yr⁻¹ in Sette et al., 2013), the amount of K that transited in the xylem sap every day was low considering the magnitude of the sap flow in the xylem for tree transpiration. This led to very low xylem sap K concentrations if a flux of K from phloem sap to xylem sap was not considered. Since values of xylem sap K concentrations measured on different plants are 1 order of magnitude higher (Nardini et al., 2010; Siebrecht et al., 2003), our results suggest that either xylem sap K concentrations are very variable between plants or an intense recirculation of K between xylem and phloem is taking place. The first hypothesis is consistent with evidence from temperate conifers that show a variation of 1 order of magnitude (the lower bound is similar to values calculated from our simulations) in terms of xylem sap K concentration during the growing season (Losso et al., 2018). However, the last hypothesis seems more plausible, since K is necessary to maintain xylem hydraulic conductivity (Oddo et al., 2011; Nardini et al., 2011) and, considering a transfer of K from phloem sap to xylem sap, the calculated K xylem sap concentrations are more in line with measurements from the literature. It also supports evidence that K⁺ ions are an essential part of many processes at the plant level (Dreyer and Michard, 2020), e.g. energy source (Dreyer et al., 2017), counter ion for NO${}_{\mathrm{3}}^{-}$, signalling (Anschütz et al., 2014), and protection against abiotic stress (Cakmak, 2005). The intense recirculation of K between xylem and phloem serves as a way to maintain homeostasis. The order of magnitude of this recirculation would have to be determined experimentally. Previous measurements have shown that up to 25 % of K⁺ ions are recirculated in tomato plants (Armstrong and Kirkby, 1979) and up to 50 % are recirculated in Ricinus communis seedlings (Marschnert et al., 1997). Results from our model suggest that the figure could be higher in eucalypt trees (maybe owing to remobilisation from wood absent in tomato plants).

5 Conclusions

The results shown here show the relevance of using a mechanistic model to explore the links between K availability and forest NPP. Our model was able to reproduce results from fertilisation experiments with a good degree of accuracy while also being able to simulate an availability gradient in K. The model allowed us to formulate the following conclusions:

1. The decrease in GPP caused by K deficiency explained most of the difference in wood productivity between K-rich and K-deficient stands. There was no need for a direct modification of the partitioning process or stoichiometric limitations to explain the observed patterns of productivity.

2. Potassium use efficiency for the production of wood increases with diminishing K availability, suggesting that the trees are able to compensate for lower K content in tissues through reduced carbon use efficiency.

3. The major importance of external inputs (weathering and deposition mainly) already observed in GPP simulations (Cornut et al., 2022) is also critical for the ecosystem's NPP.

4. The internal K fluxes in the trees that provide buffering against temporally localised K limitations are important and contribute to explaining why an adequate K supply limited to the early stage of the development of trees can be sufficient to sustain high NPPs over the entire rotation in planted forests growing in highly weathered soils.