**Development and technical paper**| 02 Jan 2023

# Accelerated photosynthesis routine in LPJmL4

Jenny Niebsch Werner von Bloh Kirsten Thonicke and Ronny Ramlau

^{1},

^{2},

^{2},

^{1}

**Jenny Niebsch et al.**Jenny Niebsch Werner von Bloh Kirsten Thonicke and Ronny Ramlau

^{1},

^{2},

^{2},

^{1}

^{1}RICAM, Altenbergerstr. 69, 4040 Linz, Austria^{2}Potsdam Institute for Climate Impact Research (PIK), Member of the Leibniz Association, 14412 Potsdam, Germany

^{1}RICAM, Altenbergerstr. 69, 4040 Linz, Austria^{2}Potsdam Institute for Climate Impact Research (PIK), Member of the Leibniz Association, 14412 Potsdam, Germany

**Correspondence**: Jenny Niebsch (jenny.niebsch@oeaw.ac.at)

**Correspondence**: Jenny Niebsch (jenny.niebsch@oeaw.ac.at)

Received: 03 May 2022 – Discussion started: 29 Jun 2022 – Revised: 14 Nov 2022 – Accepted: 21 Nov 2022 – Published: 02 Jan 2023

The increasing impacts of climate change require strategies for climate adaptation. Dynamic global vegetation models (DGVMs) are one type of multi-sectorial impact model with which the effects of multiple interacting processes in the terrestrial biosphere under climate change can be studied. The complexity of DGVMs is increasing as more and more processes, especially for plant physiology, are implemented. Therefore, there is a growing demand for increasing the computational performance of the underlying algorithms as well as ensuring their numerical accuracy. One way to approach this issue is to analyse the routines which have the potential for improved computational efficiency and/or increased accuracy when applying sophisticated mathematical methods.

In this paper, the Farquhar–Collatz photosynthesis model under water stress as implemented in the Lund–Potsdam–Jena managed Land DGVM (4.0.002) was examined.
We additionally tested the uncertainty of most important parameter of photosynthesis as an additional approach to improve model quality. We found that the numerical solution of a nonlinear equation, so far solved with the bisection method, could be significantly improved by using Newton's method instead. The latter requires the computation of the derivative of the underlying function which is presented. Model simulations show a significantly lower number of iterations to solve the equation numerically and an overall run time
reduction of the model of about 16 % depending on the chosen accuracy.
Increasing the parameters *θ* and ${\mathit{\alpha}}_{{\mathrm{C}}_{\mathrm{3}}}$ by 10 %, respectively, while keeping all other parameters at their original value, increased global gross primary production (GPP) by 2.384 and 9.542 GtC yr^{−1}, respectively. The Farquhar–Collatz photosynthesis model forms the core component in many DGVMs and land surface models. An update in the numerical solution of the nonlinear equation
in connection with adjusting globally important parameters to best known values can therefore be applied to similar photosynthesis models. Furthermore, this exercise can serve as an example for improving computationally costly routines while improving their mathematical accuracy.

- Article
(7421 KB) -
Supplement
(365 KB) - BibTeX
- EndNote

Climate change is increasingly affecting the world we live in, and that in turn affects nature's contribution to our livelihoods (Pörtner et al., 2022). Estimating the extent and impact of climate change has become more and more urgent over the last couple of decades. Earth system models (ESMs) as well as impact models are used to develop strategies for climate adaptation and mitigation to achieve the Paris climate accord (Masson-Delmotte et al., 2021; Pörtner et al., 2022). Climate change affects vegetation dynamics, biodiversity, water, and biogeochemical cycles, which could reduce the biosphere's capacity to absorb carbon from the atmosphere in the future. Dynamic global vegetation models (DGVMs) are applied to study the net effects of multiple interacting processes that affect carbon sequestration (photosynthesis) and storage (in biomass and soil), see Prentice et al. (2007). It shows the demand for reliable and consistent model projections which require continuous work on reducing model uncertainty. While increasing complexity of the models by including more and more processes in DGVMs has been matched by increasing high-performance computing capabilities over the past decades, little has been invested into identifying and optimising computationally intensive routines in the model (Reichstein et al., 2019). These routines often have a long model history as they frequently belong to the core routines stemming from the very first model version. This includes, e.g. the physiological modelling core of simulating photosynthesis in connection with atmospheric water demand or plant-water stress. The photosynthesis model is based on the Farquhar approach (Collatz et al., 1991, 1992; Farquhar et al., 1980) implemented in land surface schemes of the second generation (Pitman, 2003) followed by the first global biome models (Haxeltine and Prentice, 1996a) from which DGVMs evolved later on (Prentice et al., 2007).

The Farquhar–Collatz approach was implemented in the land surface of the SiB2 model by Sellers et al. (1992, 1996a), where it replaced their empirical photosynthesis model. The photosynthesis model in SiB2 (Sellers et al., 1996b) covers the co-limitation by Rubisco enzyme activity, light availability, and export limitation of carbon compounds. Furthermore, it covers the gradient between inner-stomatal CO_{2} concentration to the CO_{2} concentration around the leaf surface in the computation of stomatal conductance. By implementing the semi-mechanistic photosynthesis model and coupling it to transpiration via stomatal conductance, the land surface model (LSM) could then not only investigate biophysical effects of climate change but also the biogeochemical effects of rising atmospheric CO_{2} in the earth system (Pitman, 2003). The SiB2 model (Sellers et al., 1992, 1996a), the NCAR CCM2 model (Bonan et al., 1995), and the MOSES land surface model of the UK Met office (Cox et al., 1998) were among the first to implement this photosynthesis scheme and evaluate it against field campaigns. At present, the Farquhar–Collatz photosynthesis model is used in a number of land surface models of the CMIP-5 earth system models, such as the Community Atmosphere Biosphere Land Exchange (CABLE), the LSM of the Australian community climate and earth system simulator (ACCESS, see de Kauwe et al., 2015, and ref. therein), as well as the ORCHIDEE DGVM (Krinner et al., 2005) of the IPSL-CM5 earth system model (Dufresne et al., 2013). Different models of stomatal conductance were evaluated for the JSBACH LSM (Reick et al., 2013) of the Max Planck Institute earth system model (MPI-ESM) to account for hydraulic properties and drought response (Knauer et al., 2015). The Community Land Model CLM4.5
(Oleson et al., 2013) of the NCAR ESM use the Ball–Berry model of stomatal conductance and extended it to account for leaf temperature acclimation and leaf water potential (Bonan et al., 2014); a similar approach was implemented in the JULES-vn5.6 land surface model (Oliver et al., 2022) of the UK Hadley Centre ESM (Sellar et al., 2019).

While land surface models detail vertical water, energy, and carbon profiles within the canopy, which extrapolates the photosynthetic capacity calculated at the leaf level to canopy photosynthesis (Sellers et al., 1996b), stand-alone DGVMs often use a big-leaf approach and compute daytime photosynthesis for canopy conductance, which goes back to the BIOME-3 model (Haxeltine and Prentice, 1996b) that opened up the second line of vegetation models by embedding the Farquhar–Collatz photosynthesis model in a modelling framework of plant physiology and vegetation dynamics in DGVMs (Prentice et al., 2007). The Haxeltine and Prentice (1996b) implementation is used in the LPJ model family originating from Sitch et al. (2003) and the LPJ-GUESS model (Smith et al., 2001, 2014), as well as the current LPJmLv4 model (Schaphoff et al., 2018a, b). The Farquhar photosynthesis module forms the core of many other DGVMs, see e.g. Smith et al. (2001, 2014); Krinner et al. (2005). Today, 14 DGVMs (stand-alone and coupled to land surface models) contribute to the TRENDY intercomparison project (https://blogs.exeter.ac.uk/trendy/, last access: 14 December 2022), which informs the global carbon project on the state of the land carbon sink (Sitch et al., 2015).

In order to apply the model to the global land surface it is no longer sufficient to use faster or larger computing infrastructure or to try to parallelise the code as in von Bloh et al. (2010). Rather it requires the evaluation of the underlying algorithm structure of the code and, in particular, the used numerical methods. Replacing “old” numerical algorithms with modern methods will result in a significantly better run time performance while simultaneously maintaining or even increasing the accuracy of the method. We quantified the run time required by each submodule (or routine) of the LPJmL DGVM using the profiling option of the compilation command and the Linux gprof utility. We found that the repeated execution of the photosynthesis routine demands a big fraction, i.e. 38 %, of the computational time. All other routines require less than 11 %.

To illustrate our approach, our goal was to improve the computational efficiency of DGVMs by accelerating the photosynthesis module under water stress conditions using the Lund–Potsdam–Jena DGVM, Schaphoff et al. (2018a, b) as an example. A key ingredient in the modelling of photosynthesis is the determination of the ratio *λ* between intracellular and ambient CO_{2} concentration. Mathematically, *λ* is computed as a zero of a nonlinear equation *f*(*λ*)=0, which has so far been solved by a simple bisection algorithm. We expected to improve the computational efficiency by applying one of the more sophisticated solution methods, namely Regula falsi, secant and Newton's method. In this technical paper, we describe testing all three methods but found that only with Newton's method was the computational efficiency significantly improved.
Only a few detailed specialised studies mention the use of Newton's or similar methods to solve coupled balance schemes, (Collatz et al., 1991; Pearcy et al., 1997; Soo-Hyung and Lieth, 2003; Dubois et al., 2007), or extensions of the photosynthesis-transpiration scheme along the leaf–plant–soil continuum in DGVMs (Bonan et al., 2014) are mentioned, but none provide documentation on the computational efficiency or how the numerical method was implemented in the model and/or their code. In addition, we test the effect of sensitive photosynthesis parameters on the annual gross primary production (GPP) of the computationally efficient model where we build on recent work by Walker et al. (2021).

We start with a short description of the different mathematical methods to find the zeros of a general nonlinear continuous function *f* and their advantages and disadvantages. Afterwards, we introduce the relevant function *f* from the photosynthesis module and calculate its derivative. We then compare the performance of Newton's algorithm and bisection in terms of the number of iterations and the computational time that is necessary to achieve a given accuracy. Finally, we benchmark the updated LPJmL version to show that the simulated vegetation dynamics as well as storage and fluxes of carbon and water remain robust.

The computation of the ratio *λ* between intracellular and ambient CO_{2} concentrations requires us to compute the zero of a function *f*(*λ*). In most cases, this task cannot be solved analytically but requires a numerical approach, mostly based on iterative methods.
Given a nonlinear continuous function $f:\mathbb{R}\to \mathbb{R}$, we want to find the zero(s) *x*_{s} of this function within a certain interval [*a*,*b*]. While bisection, regula falsi, and secant methods are very simple to implement, Newton's method requires the computation of the derivative of *f*, which will be provided for the photosynthesis equation described in Sect. 3.2.

Here, the computational efficiency is determined by the speed of convergence. To compare the methods with respect to the speed of convergence we define the order of convergence as follows:
let *x*_{s} be a zero of *f* found by computing a sequence (*x*_{k}) of approximate solutions via an iteration scheme. The iteration method has the order of convergence *p* if

with $\mathrm{0}<K<\mathrm{\infty}$ and *K*<1 for *p*=1.
Thus a high order of convergence implies a fast convergence, which on the other hand means fewer iteration steps.
Numerically, the iteration is stopped either if the function value *f*(*x*_{k}) of the iterate *x*_{k} is almost zero, i.e. less than a given accuracy *y*_{acc}, or if the iterate itself changes less than a given accuracy $|{x}_{k}-{x}_{k-\mathrm{1}}|<{x}_{\mathrm{acc}}$.

Let us introduce some of the methods in the following subsections, see Schwarz and Köckler (2009) for details.

## 2.1 Bisection

For bisection we have to choose [*a*,*b*] such that $f\left(a\right)\cdot f\left(b\right)<\mathrm{0}$, i.e. *f*(*a*) and *f*(*b*) have different signs.
We compute the midpoint of the interval ${x}_{\mathrm{m}}=\frac{a+b}{\mathrm{2}}$ and its function value *f*(*x*_{m}). If $\left|f\right({x}_{\mathrm{m}}\left)\right|<{y}_{\mathrm{acc}}$ the search is complete, if not we check if $f\left(a\right)\cdot f\left({x}_{\mathrm{m}}\right)<\mathrm{0}$. If the latter is the case, *x*_{s} has to be in the interval [*a*,*x*_{m}] or otherwise in [*x*_{m},*b*]. We repeat this bisection until either $\left|f\right({x}_{k}\left)\right|<{y}_{\mathrm{acc}}$ or $|{x}_{k}-{x}_{k-\mathrm{1}}|<{x}_{\mathrm{acc}}$. This method always converges but slowly with convergence order *p*=1, i.e. linear convergence.

## 2.2 Regula falsi

For the regula falsi method, we also need to choose *a*,*b* such that $f\left(a\right)\cdot f\left(b\right)<\mathrm{0}$. Instead of the midpoint of [*a*,*b*], we compute the next iterate *x*_{1} for an approximation of *x*_{s} by computing the zero of the linear function through the points (*a*|*f*(*a*)) and (*b*|*f*(*b*)). Again we check if $\left|f\right({x}_{\mathrm{1}}\left)\right|<{y}_{\mathrm{acc}}$ and abort or check if $f\left(a\right)\cdot f\left({x}_{\mathrm{1}}\right)<\mathrm{0}$, and repeat this procedure either with [*a*,*x*_{1}] or [*x*_{1},*b*]. Convergence is always assured and is also linear, i.e. *p*=1.

## 2.3 Secant method

The secant method only differs from the regula falsi in that the starting values *a*=*x*_{0} and *b*=*x*_{1} do not have to fulfil the condition $f\left(a\right)\cdot f\left(b\right)<\mathrm{0}$.
The next iterate is computed by

This method can fail to converge depending on the starting values. If the method converges, it does so with order *p*=1.618. Since the conditions on the starting values to ensure convergence depend on the knowledge of *x*_{s}, in practise *a* and *b* still have to fulfil the condition $f\left(a\right)\cdot f\left(b\right)<\mathrm{0}$.

## 2.4 Newton's method

Newton's method starts at an arbitrary approximation *x*_{0} of *x*_{s} and uses the tangent of the function *f* at (*x*_{0},*f*(*x*_{0})) to compute the next iterate *x*_{1} as the zero of the tangent. This is repeated, thus the next iterate is always computed from the previous one by

provided that ${f}^{\prime}\left({x}_{k}\right)\ne \mathrm{0}$. The method belongs to the class of fixed point iterations because the computation of the next iterate depends on the previous iterate only. If *f* is three times differentiable on [*a*,*b*] and ${f}^{\prime}\left({x}_{s}\right)\ne \mathrm{0}$, then there exists an interval $I=[{x}_{s}-\mathit{\delta},{x}_{s}+\mathit{\delta}]$ such that *f* is a contraction on *I*. It implies that for every start value from *I*, the method converges at least with order *p*=2 (Schwarz and Köckler, 2009).
We remark that the gain in convergence speed has to be weighted against the time it takes to compute the derivative of *f*.

We now analyse the difference in speed of convergence between the bisection and Newton's methods when applied to the optimisation equation of the photosynthesis routine of the LPJmL DGVM.

## 3.1 Definition of the function *f*

In presenting the function *f*(*λ*), we follow the nomenclature of Schaphoff et al. (2018a), which contains a detailed description of the derivation of this function. A list of the used symbols is given in Appendix A.
We want to find $\mathit{\lambda}=\frac{{c}_{i}}{{c}_{\mathrm{a}}}=\frac{{p}_{i}}{{p}_{\mathrm{a}}}$, i.e. the ratio between the intracellular and ambient CO_{2} concentration, or partial pressure, respectively, as the solution of the following equation:

Here *A*_{nd} is the net daily photosynthesis, *R*_{leaf} is the leaf respiration, dayl is the hours of daylight, *p*_{a} is the ambient partial pressure, *g*_{c} is the canopy conductance, and *g*_{min} is the minimum canopy conductance for a specific plant functional type (PFT).
The first term is the photosynthesis during daylight.
It is the gross daily photosynthesis *A*_{gd} minus leaf respiration, ${A}_{\mathrm{nd}}\left(\mathit{\lambda}\right)={A}_{\mathrm{gd}}\left(\mathit{\lambda}\right)-{R}_{\mathrm{leaf}}$. The second term represents the dark respiration, i.e. respiration during night time. The third term represents the photosynthesis that is possible to achieve a potential canopy conductance. In finding *λ* such that *f*(*λ*)≈0 we actually balance both light- and Rubisco-limited photosynthesis (first two terms) and photosynthesis related to the potential canopy conductance.

To shorten the formulas we define the abbreviation ${C}_{\mathrm{pg}}=\frac{{p}_{\mathrm{a}}({g}_{\mathrm{c}}-{g}_{\mathrm{min}})}{\mathrm{1.6}}$ as

The second summand does not depend on *λ*, whereas *A*_{gd}(*λ*) has a more complex representation. The gross photosynthesis rate *A*_{g} is the minimum of the light-limited, *J*_{C}, and Rubisco-limited photosynthesis rate, *J*_{E}. It can be shown that the minimum can be computed as

where *θ* is a shape parameter that allows for a gradual transition from one limitation to the other.

Light-limited photosynthesis depends on the absorbed photosynthetically active radiation (APAR); Rubisco-limited photosynthesis is determined by the maximum Rubisco capacity *V*_{m}:

Setting the internal partial pressure
*p*_{i}=*λ**p*_{a} and using another abbreviation ${C}_{K}={K}_{\mathrm{c}}(\mathrm{1}+\frac{\left[{\mathrm{O}}_{\mathrm{2}}\right]}{{K}_{\mathrm{O}}})$, where *K*_{c} is the Michaelis constant for CO_{2} and [O_{2}] and *K*_{O} are the partial pressure and the Michaelis constant for oxygen, we have

Here, ${\mathit{\alpha}}_{{\mathrm{C}}_{\mathrm{3}}}$ and ${\mathit{\alpha}}_{{\mathrm{C}}_{\mathrm{4}}}$ are the intrinsic quantum efficiencies for CO_{2} uptake in C_{3} and C_{4} plants, respectively. Γ_{*} is the carbon dioxide compensation point and *T*_{stress} is a temperature stress function defined as

with *T*_{d} as the daily air temperature and *T*_{1} to *T*_{4} being PFT-specific temperature parameters (Sitch et al., 2000). LPJmL simulates vegetation dynamics for the 10 PFTs; we provide the parameter values used for *T*_{1} to *T*_{4} in Appendix A, Table A1, for the PFT types from Schaphoff et al. (2018a).

## 3.2 Derivative of *f*

To compute the derivative *f*^{′} of *f* we rearrange Eq. (5):

Since the last two terms are constant the derivative is given by

To determine ${A}_{\mathrm{gd}}^{\prime}$ we apply sum, chain, and product rule of differentiation to Eq. (6) and get

The derivatives of *J*_{E} and *J*_{C} are given by

To compute ${C}_{\mathrm{1}}^{\prime}$ from Eq. (9) and ${C}_{\mathrm{2}}^{\prime}$ from Eq. (10) we use the quotient rule

We describe the consequent changes in the model code which were required to implement the computation of the derivative fcnd(*λ*) in the Appendix B.

The function *f* is defined for all *λ*>0, as long as $\left({J}_{E}\right(\mathit{\lambda})+{J}_{C}(\mathit{\lambda}){)}^{\mathrm{2}}\ge \mathrm{4}\mathit{\theta}{J}_{E}(\mathit{\lambda}\left){J}_{C}\right(\mathit{\lambda})$. As a composition of at least 3 times differentiable functions, it fulfils the differentiability condition of Newton's method. The parameters in the definition of *f* vary with the geographic location and season. A plot of *f* for parameters from different locations (boreal, temperate, and tropical) and at different times can be seen in Fig. 1.

The condition ${f}^{\prime}\left(\mathit{\lambda}\right)\ne \mathrm{0}$ as well as the suitability of a staring value can not be generally ensured. In all our computations convergence was not a problem. To be on the safe side, one can implement a hybrid method that switches to bisection if convergence of the iterates does not occur.

We have tested the different methods in the routine regarding computational time and number of iterations for given accuracy *x*_{acc}. There was no significant speed up with the secant and regula falsi method. Hence, we concentrated on the comparison of bisection and Newton's methods and describe the outcome in this section.

In a first test, the LPJmL model was run over 120 simulation years and the number of iterations in the bisection and Newton's routine was counted and averaged over all grid cells and one year (Fig. 2). For *x*_{acc}=0.01 this number was about 3 for Newton's method and 7 for bisection (dotted lines in Fig. 2). When *x*_{acc} was set to 0.001 the number of iterations with Newton's method increased only slightly, whereas the bisection method needed 9 to 10 iterations (solid lines in Fig. 2). Until now, the bisection algorithm used 10 as the maximal number of iterations. Using maximum 10 iterations fits into the interval width of ${\mathrm{2}}^{-\mathrm{10}}\approx \mathrm{0.001}$, our accuracy measure *x*_{acc}. Increasing the maximum number of iterations had no effect on the number of required iterations. We conclude that Newton's method reduces the necessary number of iteration to a third.

In a next step, a spin-up run of LPJmL over 5000 simulation years was conducted to compare the time performance using both routines. Usually, LPJmL simulation experiments start from bare ground, i.e. initial vegetation conditions are not prescribed. Therefore, a spin-up run is used to bring all vegetation and soil carbon pools into equilibrium with climate.
For the usually implemented accuracy *x*_{acc}=0.1 the computation time for 5000 years was about 5250 s in both cases. This means that the advantage of Newton's method in terms of iteration numbers is levelled by the additional time for computing the derivative of *f*. For *x*_{acc}=0.01, the bisection method needed 6700 s, while Newton's method needed 5600 s. Thus a reduction of about 16 % in time could be observed. It implies that with almost the same amount of time (5250 s vs. 5600 s) a higher accuracy can be achieved with Newton's method (Fig. 3). While the accuracy *y*_{acc} does not increase significantly for the bisection method for *x*_{acc}=0.001, we gain a 2 orders of magnitude increase in *y*_{acc} for the Newton's method. As a result, a change of *x*_{acc} from 0.1 to 0.01 will be permanently implemented in the LPJmL model for future model applications. We expect that with the implementation of new model developments that affect the photosynthesis module (e.g. nutrient limitation from nitrogen and leaf temperatures) an efficient and increased model accuracy (*y*_{acc}) for finding the zero of *f*(*λ*) will be even more important. It can be expected that the computation time for the bisection method would increase substantially, while increasing only moderately for Newton's method.

In order to check if the implementation of Newton's method is robust for all important model variables, we performed a transient simulation with the LPJmL model starting from the spin-up and covering the years 1901–2000. Model configuration and input data are as in Schaphoff et al. (2018a). We compared the main diagnostic variables of the published LPJmL4.0 version against the version using Newton's method (see Appendix C). We found that most global diagnostic variables related to fluxes and storage of carbon and water had differences of $<\pm \mathrm{1.0}$ %, including total vegetated area. Only marginal changes (+3 gC per m^{2} and month) in net primary productivity (NPP), heterotrophic respiration, and evaporation are seen mainly in Europe and southern as well as southeastern Asia. The reductions in carbon storage in litter and soil are very small and apply only to the boreal zone across the Northern Hemisphere and central Europe (compare spatial maps of carbon and water variables in Appendix C).

The photosynthesis module is also applied to the crop functional types and managed grassland within LPJmL4.0. Therefore, sawing dates, crop productivity, and harvest are among the simulated variables. Comparing both model versions in the model benchmark, we found that global harvest changed for a number of crops. Rainfed and irrigated rice increased by 5 % and 8 %, respectively, mainly in India and southeast Asia. Harvest of rainfed temperate cereals increased by <1.0 %, mainly found in central Europe. Harvest of irrigated temperate cereals (incl. wheat) increased by 4.5 %, which mainly applied to India as well. Harvest of irrigated and rainfed soybeans increased by 2.3 % and 1.5 % globally; the differences are mainly found in the US and Brazil. All other crop functional types had marginal to zero changes in global productivity as well as simulated harvest (see Table in Appendix C).

For all global carbon pools (vegetation and soil) and carbon (GPP, heterotrophic respiration, and fire emissions) as well as water fluxes (transpiration and runoff) we found no difference in the temporal changes in the transient simulation over the 20th century. All variables showed similar, if not identical, dynamics (data not shown).
Small changes were found in the fractional coverage of plant functional types, i.e. most differences were negligible. The fractional coverage of temperate broadleaved summergreen trees increased by 4.8 % globally, which mainly applies to Europe, the northeastern USA, and parts of China. Increases in temperate C_{3} grasses are found in the boreal zone, summing up to 4.8 % globally. Marginal changes of <0.5 *%* per grid cell are found for all other PFTs, which imply small adjustments in vegetation composition in these vegetation zones (see difference maps in Appendix C).
Comparisons using flux tower measurements on carbon and water fluxes as well as discharge data showed no differences so we can conclude that also for these variables the results are robust (data not shown). We can therefore conclude that the LPJmL results were robust before but are now achieved due to improved accuracy of the photosynthesis routine.

After improving the computational efficiency and numerical precision, we can now test the parameter uncertainties following Walker et al. (2021), who tested the sensitivity of *θ*, ${\mathit{\alpha}}_{{\mathrm{C}}_{\mathrm{3}}}$, ${b}_{{\mathrm{C}}_{\mathrm{3}}}$, *k*_{c25}, and *K*_{o25} on their impacts on global GPP. The LPJmL model computes *V*_{m} as follows (Schaphoff et al., 2018a, Eq. 35):

Therefore, the sensitivity of *V*_{cmax} results from varying ${b}_{{\mathrm{C}}_{\mathrm{3}}}$ indirectly since the reciprocal of ${b}_{{\mathrm{C}}_{\mathrm{3}}}$ is used to calculate *V*_{cmax} in a linear equation. Varying ${b}_{{\mathrm{C}}_{\mathrm{3}}}$ is therefore the adequate sensitivity test which relates to *V*_{cmax}. We varied each parameter by 10 % independently and find that *θ* (${\mathit{\alpha}}_{{\mathrm{C}}_{\mathrm{3}}}$, ${b}_{{\mathrm{C}}_{\mathrm{3}}}$, *k*_{c25}, *K*_{o25}) increases global annual GPP (AGPP, hereafter) by 1.67 % (+6.69 %, −1.67 %, −0.35 %, and +0.14 %). Table 1 shows the difference of the two most important parameter on global AGPP.

Geographically, increasing *θ* yields higher AGPP mainly in the tropics and temperate forest regions, where AGPP increases up to 100 gC m^{−2}. However, AGPP increases between 200 and 500 gC m^{−2} when changing ${\mathit{\alpha}}_{{\mathrm{C}}_{\mathrm{3}}}$, see Fig. 4. It turns out that AGPP is increased in all regions, where LPJmL simulates woody PFTs. Also here, the largest effects are seen in (sub-)tropical and temperate regions which span larger areas than the areas with increased AGPP as a result of varying *θ*.

We remark that future work on the photosynthesis approach could focus on the new Johnson and Berry scheme (Johnson and Berry , 2021) with the advantage of calculating gas exchange and relying less on empirical coefficients.

The computational load of dynamic global vegetation models, caused by increased complexity of the modelling processes, has so far been counteracted by the high-performance computing systems used. However, more recently it has become clear that updates in computing infrastructure are not sufficient anymore. Consequently, we proposed to carefully evaluate the algorithmic structure of DGVMs and identify and update routines that can benefit from the use of modern mathematical methods. As a showcase, we investigated the photosynthesis model in the LPJmL DGVM. Specifically, we investigated the computation of the ratio *λ* between intracellular and ambient CO_{2}, which is obtained as the zero of a function *f*. We proposed to replace the so far used bisection method with a Newton method, which is known to converge significantly faster. We carefully compared the model performance of the published LPJmL4.0 version with the version developed in this study and found that the model performance is robust. Using a more sophisticated mathematical method in the photosynthesis module allowed for a higher precision in the computation of *λ* and resulted in slightly increased productivity in continental and mountainous areas. We think that the new results are more accurate than the previous version due to the higher accuracy of the Newton method visible in Fig. 3. With the currently implemented accuracy bounds, the run time of the model with the Newton routine implemented is about 16 % lower than the old version. This advantage will be much more prominent if the complexity of the model is further extended or if more accurate modelling results are required. Consequently, the Newton-based routine will be implemented in the LPJmL model. Additionally, we believe that the Newton method can also be applied to photosynthesis modules in other DGVMs and can increase model accuracy and/or computational efficiency.

General parameters used in the photosynthesis routine. PFT is plant functional type.

A_{nd} |
daily net photosynthesis |

dayl | day length |

R_{leaf} |
leaf respiration |

p_{a} |
ambient partial pressure |

g_{c} |
canopy conductance |

g_{min} |
PFT-specific minimum canopy conductance |

A_{gd} |
daily gross photosynthesis |

θ |
co-limitation (shape) parameter |

J_{E} |
light-limited photosynthesis rate |

J_{C} |
Rubisco-limited photosynthesis rate |

APAR | absorbed photosynthetically active radiation |

V_{m} |
maximum Rubisco capacity |

K_{C} |
Michaelis constant for CO_{2} |

[O_{2}] |
O_{2} partial pressure |

K_{O} |
Michaelis constant for O_{2} |

T_{stress} |
temperature stress function limiting photosynthesis |

at low and high temperatures | |

${\mathit{\alpha}}_{{\mathrm{C}}_{\mathrm{3}}}$ | intrinsic quantum efficiencies for CO_{2} uptake in C_{3} plants |

${\mathit{\alpha}}_{{\mathrm{C}}_{\mathrm{4}}}$ | intrinsic quantum efficiencies for CO_{2} uptake in C_{4} plants |

Γ_{*} |
carbon dioxide compensation point |

${\mathit{\lambda}}_{{\mathrm{maxC}}_{\mathrm{4}}}$ | maximum ratio of intracellular to ambient CO_{2} for C_{4}-photosynthesis |

To implement Newton's method in the LPJmL code, changes had to be made in the functions `photosynthesis.c`

, `gp_sum.c`

, and ` water_stressed.c`

. (separate file)

New function `newton.c`

: see source code in a separate file.

*Remark*.
The function `photosynthesis.c`

within LPJmL computes the value ${A}_{\mathrm{nd}}\left(\mathit{\lambda}\right)+\left(\mathrm{1}-\frac{\mathrm{dayl}}{\mathrm{24}}\right){R}_{\mathrm{leaf}}$ for a given *λ*.
In the function `water_stressed.c`

the function `fcn`

(*λ*) is defined as
$\mathrm{fcn}\left(\mathit{\lambda}\right)={C}_{\mathrm{pg}}\times (\mathrm{1}-\mathit{\lambda})-\mathrm{photosythesis}\left(\mathit{\lambda}\right)$, i.e. $\mathrm{fcn}=-f$.
In order to use Newton's method we have to compute not only fcn(*λ*) but also its derivative $\mathrm{fcnd}\left(\mathit{\lambda}\right)=-{f}^{\prime}\left(\mathit{\lambda}\right)$.

The benchmark table of global status variables (Table C1) compares two model versions against each other and to literature values were available. The following Figs. D1–D6 show globally important variables simulated using the Newton approach (benchmark run) and the bisection method (run) as time series and maps.

Literature: ^{a} Olson et al. (1985). ^{b} Saugier et al. (2001). ^{c} WBGU (1998). ^{d} Batjes (1996). ^{e} Eswaran et al. (1993). ^{f} Post et al. (1982). ^{g} Seiler and Crutzen (1980). ^{h} Andreae and Merlet (2001). ^{i} Ito and Penner (2004). ^{j} van der Werf et al. (2004). ^{k} Vitousek et al. (1986). ^{l} Ramakrishna et al. (2003). ^{m} FAOSTAT (2009).

The model code is available at https://doi.org/10.5281/zenodo.6644541 (Niebsch et al., 2022).

The supplement related to this article is available online at: https://doi.org/10.5194/gmd-16-17-2023-supplement.

JN and RR performed the mathematical analysis, JN and WvB implemented and tested the new numerical methods, and WvB conducted the simulation experiments and analysed the model performance and computation efficiency. JN and KT wrote the paper and all authors contributed to the writing of the paper and discussion of the model study throughout to develop the work.

The contact author has declared that none of the authors has any competing interests.

Publisher’s note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

The authors gratefully acknowledge the European Regional Development Fund (ERDF), the German Federal Ministry of Education and Research, and the Land Brandenburg for supporting this project by providing resources on the high-performance computer system at the Potsdam Institute for Climate Impact Research. We thank Marie Hemmen from PIK for her support in benchmarking the LPJmL model.

The high-performance computing system at PIK was funded by the European Regional Development Fund (ERDF), the German Federal Ministry of Education and Research, and the Land Brandenburg.

This paper was edited by Carlos Sierra and reviewed by two anonymous referees.

Andreae, M. O. and Merlet, P.: Emission of trace gases and aerosols from biomass burning, Global Biogeochem. Cy., 15, 955–966, https://doi.org/10.1029/2000GB001382, 2001. a

Batjes, N. H.: Total carbon and nitrogen in the soils of the world, Eur. J. Soil Sci., 47, 151–163, https://doi.org/10.1111/j.1365-2389.1996.tb01386.x, 1996. a

Bonan, G. B.: Land–atmosphere CO_{2} exchange simulated by a land surface process model coupled to an atmospheric general circulation model, J. Geophys. Res., 100, 2817–2831, 1995. a

Bonan, G. B., Williams, M., Fisher, R. A., and Oleson, K. W.: Modeling stomatal conductance in the earth system: linking leaf water-use efficiency and water transport along the soil–plant–atmosphere continuum, Geosci. Model Dev., 7, 2193–2222, https://doi.org/10.5194/gmd-7-2193-2014, 2014. a, b

Collatz, G. J., Ball, J. T., Grivet, C., and Berry, J. A.: Physiological and environmental regulation of stomatal conductance, photosynthesis and transpiration: a model that includes a laminar boundary layer, Agric. For. Meteorol., 54, 107–136, https://doi.org/10.1016/0168-1923(91)90002-8, 1991. a, b

Collatz, G. J., Ribas-Carbo, M., and Berry, J. A.: Coupled Photosynthesis-Stomatal Conductance Model for Leaves of C_{4} Plants, Funct. Plant Biol., 19, 519–538, https://doi.org/10.1071/PP9920519, 1992. a

Cox, P. M., Huntingford, C., and Harding, R. J.: A canopy conductance and photosynthesis model for use in a GCM land surface scheme, J. Hydrol., 212–213, 79–94, 1998. a

De Kauwe, M. G., Kala, J., Lin, Y.-S., Pitman, A. J., Medlyn, B. E., Duursma, R. A., Abramowitz, G., Wang, Y.-P., and Miralles, D. G.: A test of an optimal stomatal conductance scheme within the CABLE land surface model, Geosci. Model Dev., 8, 431–452, https://doi.org/10.5194/gmd-8-431-2015, 2015. a

Dubois, J.-J. B., Fiscus, E. L., Booker, F. L., Flowers, M. D., and Reid, C. D.: Optimizing the statistical estimation of the parameters of the Farquhar–von Caemmerer–Berry model of photosynthesis, New Phytol., 176, 402–414, https://doi.org/10.1111/j.1469-8137.2007.02182.x, 2007. a

Dufresne, J. L., Foujols, M. A., and Denvil, S., et al.: Climate change projections using the IPSL-CM5 Earth System Model: from CMIP3 to CMIP5, Clim. Dynam., 40, 2123–2165, https://doi.org/10.1007/s00382-012-1636-1, 2013. a

Eswaran, H., Van Den Berg, E., and Reich, P.: Organic Carbon in Soils of the World, Soil Sci. Soc. Am. J., 57, 192–194, https://doi.org/10.2136/sssaj1993.03615995005700010034x, 1993. a

FAOSTAT: FAOSTAT database, Food and Agriculture Organization of the United Nations, Rome, http://www.fao.org/faostat/en/ (last access: 28 February 2009), 2009. a

Farquhar, G. D., von Caemmerer, S., and Berry, J. A.: A biochemical
model of photosynthetic CO_{2} assimilation in leaves of C_{3} species,
Planta, 149, 78–90, https://doi.org/10.1007/BF00386231, 1980. a

Haxeltine, A. and Prentice, I. C.: BIOME3: An equilibrium terrestrial biosphere model based on ecophysiological constraints, resource availability, and competition among plant functional types, Global Biogeochem. Cy., 10, 693–709, https://doi.org/10.1029/96GB02344, 1996a. a

Haxeltine, A. and Prentice, I. C.: A general model for the light-use efficiency of primary production, Funct. Ecol., 10, 551–561, https://doi.org/10.2307/2390165, 1996b. a, b

Ito, A. and Penner, J. E.: Global estimates of biomass burning emissions based on satellite imagery for the year 2000, J. Geophys. Res.-Atmos., 109, D14S05, https://doi.org/10.1029/2003JD004423, 2004. a

Johnson, J. E. and Berry, J. A.: The role of cytochrome b6f in the control of steady- state photosynthesis: a conceptual and quantitative model, Photosynth. Res., 148, 101–136, 2021. a

Knauer, J., Werner, C., and Zaehle, S.: Evaluating stomatal models and their atmospheric drought response in a land surface scheme: A multibiome analysis, J. Geophys. Res.-Biogeo., 120, 1894–1911, https://doi.org/10.1002/2015JG003114, 2015. a

Krinner, G., Viovy, N., Noblet-Ducoudré, N. D., Ogée, J., Polcher, J., Friedlingstein, P., Ciais, P., Sitch, S., and Prentice, I. C.: A dynamic global vegetation model for studies of the coupled atmosphere-biosphere system, Global Biogeochem. Cy., 19, GB1015, https://doi.org/10.1029/2003GB002199, 2005. a, b

Masson-Delmotte, V., Zhai, P., Pirani, A., Connors, S.L., Péan, C., Berger, S., Caud, N., Chen, Y., Goldfarb, L., Gomis, M.I., Huang, M., Leitzell, K., Lonnoy, E., Matthews, J. B. R., Maycock, T. K., Waterfield, T., Yelekçi, O., Yu, R., and Zhou, B. (Eds.): IPCC, 2021: Summary for Policymakers, in: Climate Change 2021: The Physical Science Basis. Contribution of Working Group I to the Sixth Assessment Report of the Intergovernmental Panel on Climate Change, Cambridge University Press, https://doi.org/10.1017/9781009157896, 2021. a

Niebsch, J., von Bloh, W., Thonicke, K., and Ramlau, R.: LPJmL Version 4 with Newton root finding method (4.00.4), Zenodo [code and data set], https://doi.org/10.5281/zenodo.6644541, 2022. a

Oleson, K. W., Lawrence, D. M., Bonan, G. B., Drewniak, B., Huang, M., Koven, C. D., Levis, S., Li, F., Riley, W. J., Subin, Z. M., Swenson, S. C., Thornton, P. E., Bozbiyik, A., Fisher, R., Heald, C. L., Kluzek, E., Lamarque, J.-F., Lawrence, P. J., Leung, L. R., Lipscomb, W., Muszala, S., Ricciuto, D. M., Sacks, W., Sun, Y., Tang, J., and Yang, Z.-L.: Technical description of version 4.5 of the Community Land Model (CLM), NCAR Tech. Note NCAR/TN-503+STR, National Center for Atmospheric Research, Boulder, Colorado, 420 pp., 2013. a

Olson, J. S., Watts, J. A., and Allison, L. J.: Major world ecosystem complexes ranked by carbon in live vegetation: A Database, NDP-017. Carbon Dioxide Information Center, edited by: O. R. N. L., Oak Ridge, Tennessee, https://doi.org/10.3334/CDIAC/lue.ndp017, 1985. a

Oliver, R. J., Mercado, L. M., Clark, D. B., Huntingford, C., Taylor, C. M., Vidale, P. L., McGuire, P. C., Todt, M., Folwell, S., Shamsudheen Semeena, V., and Medlyn, B. E.: Improved representation of plant physiology in the JULES-vn5.6 land surface model: photosynthesis, stomatal conductance and thermal acclimation, Geosci. Model Dev., 15, 5567–5592, https://doi.org/10.5194/gmd-15-5567-2022, 2022. a

Pearcy, R. W., Gross, L. J., and He, D.: An improved dynamic model of photosynthesis for estimation of carbon gain in sunfleck light regimes, Plant Cell Environ., 20, 411–424, https://doi.org/10.1046/j.1365-3040.1997.d01-88.x, 1997. a

Pitman, A. J.: The evolution of, and revolution in, land surface schemes designed for climate models, Int. J. Climatol., 23, 479–510, https://doi.org/10.1002/joc.893, 2003. a, b

Pörtner, H.-O., Roberts, D. C., Poloczanska, E. S., Mintenbeck, K., Tignor, M., Alegria, A., Craig, M., Langsdorf, S., Löschke, S., Möller, V., Okem, A. (Eds.): IPCC, 2022: Summary for Policymakers In: Climate Change 2022: Impacts, Adaptation, and Vulnerability. Contribution of Working Group II to the Sixth Assessment Report of the Intergovernmental Panel on Climate Change, Cambridge University Press, https://doi.org/10.1017/9781009325844, 2022. a, b

Post, W. M., Emanuel, W. R., Zinke, P. J., and Stangenberger, A. G.: Soil carbon pools and world life zones, Nature, 298, 156–159, https://doi.org/10.1038/298156a0, 1982. a

Prentice, I. C., Bondeau, A., Cramer, W., Harrison, S. P., Hickler, T., Lucht, W., Sitch, S., Smith, B., and Sykes, M. T.: Dynamic global vegetation modelling: quantifying terrestrial ecosystem responses to large-scale environmental change, in: Terrestrial Ecosystems in a Changing World, edited by: Canadell, J. G., Pataki, D. E., and Pitelka L. F., Springer, Springer Nature, https://doi.org/10.1007/978-3-540-32730-1, 2007. a, b, c

Ramakrishna R. Nemani, Keeling, C. D., Hashimoto, H., Jolly, W. M., Piper, S. C., Tucker, C. J., Myneni, R. B., and Running, S. W.: Climate-Driven Increases in Global Terrestrial Net Primary Production from 1982 to 1999, Science, 300, 1560–1563, https://doi.org/10.1126/science.1082750, 2003. a

Reichstein, M., Camps-Valls, G., Stevens, B., Jung, M., Denzler, J., Carvalhais, N., and Prabhat: Deep learning and process understanding for data-driven Earth system science, Nature, 566, 195–204, https://doi.org/10.1038/s41586-019-0912-1, 2019. a

Reick, C., Raddatz, T., Brovkin, V., and Gayler, V.: Representation of natural and anthropogenic land cover change in MPI-ESM, J. Adv. Model. Earth Sy., 5, 459–482, https://doi.org/10.1002/jame.20022, 2013. a

Saugier, B., Roy, J., and Mooney, H. A.: Estimations of Global Terrestrial Productivity: Converging toward a Single Number?, in: Terrestrial Global Productivity, edited by: Roy, J., Saugier, B., and Mooney, H. A., Academic Press, San Diego, 543–555, ISBN: 0-12-505290-1, 2001. a

Schaphoff, S., Forkel, M., Müller, C., Knauer, J., von Bloh, W., Gerten, D., Jägermeyr, J., Lucht, W., Rammig, A., Thonicke, K., and Waha, K.: LPJmL4 – a dynamic global vegetation model with managed land – Part 2: Model evaluation, Geosci. Model Dev., 11, 1377–1403, https://doi.org/10.5194/gmd-11-1377-2018, 2018a. a, b, c, d, e, f, g

Schaphoff, S., Forkel, M., Müller, C., Knauer, J., von Bloh, W., Gerten, D., Jägermeyr, J., Lucht, W., Rammig, A., Thonicke, K., and Waha, K.: LPJmL4 – a dynamic global vegetation model with managed land – Part 2: Model evaluation, Geosci. Model Dev., 11, 1377–1403, https://doi.org/10.5194/gmd-11-1377-2018, 2018b. a, b

Schwarz, H. R. and Köckler, N.: Numerische Mathematik, 7th edn., Vieweg + Teubner, Wiesbaden, https://doi.org/10.1007/978-3-8348-9282-9, 2009. a, b

Seiler, W. and Crutzen, P. J.: Estimates of gross and net fluxes of carbon between the biosphere and the atmosphere from biomass burning, Climatic Change, 2, 207–247, https://doi.org/10.1007/BF00137988, 1980. a

Sellar, A. A., Jones, C. G., Mulcahy, J. P., Tang, Y., Yool, A., Wiltshire, A., O'Connor, F. M., Stringer, M., Hill, R., Palmieri, J., Woodward, S., de Mora, L., Kuhlbrodt, T., Rumbold, S. T., Kelley, D. I., Ellis, R., Johnson, C. E., Walton, J., Abraham, N. L., Andrews, M. B., Andrews, T., Archibald, A. T., Berthou, S., Burke, E., Blockley, E., Carslaw, K., Dalvi, M., Edwards, J., Folberth, G. A., Gedney, N., Griffiths, P. T., Harper, A. B., Hendry, M. A., Hewitt, A. J., Johnson, B., Jones, A., Jones, C. D., Keeble, J., Liddicoat, S., Morgenstern, O., Parker, R. J., Predoi, V., Robertson, E., Siahaan, A., Smith, R. S., Swaminathan, R., Woodhouse, M. T., Zeng, G., and Zerroukat, M.: UKESM1: Description and Evaluation of the U.K. Earth System Model, J. Adv. Model. Earth Sy., 11, 4513–4558, https://doi.org/10.1029/2019MS001739, 2019. a

Sellers, P. J., Berry, J. A., Collatz, G. J., Field, C. B., and Hall, F. G.: Canopy reflectance, photosynthesis, and transpiration. III. A reanalysis using improved leaf models and a new canopy integration scheme, Remote Sens. Environ., 42, 187–216, https://doi.org/10.1016/0034-4257(92)90102-P, 1992. a, b

Sellers, P. J., Randall, D. A., Collatz, G. J., Berry, J. A., Field, C. B., Dazlich, D. A., Zhang, C., Collelo, G. D., and Bounoua, L.: A revised land surface parameterization (SiB2) for atmospheric GCMs. Part I: Model formulation, J. Climate, 9, 676–705, https://doi.org/10.1175/1520-0442(1996)009<0676:ARLSPF>2.0.CO;2, 1996a. a, b

Sellers, P. J., Tucker, C. J., Collatz, G. J., Los, S. O., Justice, C. O., Dazlich, D. A., and Randall, D. A.: A revised land surface parameterization (SiB2) for atmospheric GCMs. Part II: The generation of global fields of terrestrial biophysical parameters from satellite data, J. Climate, 9, 706–737, https://doi.org/10.1175/1520-0442(1996)009<0706:ARLSPF>2.0.CO;2, 1996b. a, b

Sitch, S., Prentice, I. C., Smith, B., Cramer, W., Kaplan, J., Lucht, W., Sykes, M., Thonicke, K., and Venevsky, S.: LPJ- a coupled model of vegetation dynamics and the terrestrial carbon cycle, Doctoral dissertation, Institute of Plant Ecology, Lund University, Lund, 213 pp., 2000. a

Sitch, S., Smith, B., Prentice, I. C., Arneth, A., Bondeau, A., Cramer, W., Kaplan, J. O., Levis, S., Lucht, W., Sykes, M. T., Thonicke, K., and Venevsky, S.: Evaluation of ecosystem dynamics, plant geography and terrestrial carbon cycling in the LPJ dynamic global vegetation model, Glob. Change Biol., 9, 161–185, https://doi.org/10.1046/j.1365-2486.2003.00569.x, 2003. a

Sitch, S., Friedlingstein, P., Gruber, N., Jones, S. D., Murray-Tortarolo, G., Ahlström, A., Doney, S. C., Graven, H., Heinze, C., Huntingford, C., Levis, S., Levy, P. E., Lomas, M., Poulter, B., Viovy, N., Zaehle, S., Zeng, N., Arneth, A., Bonan, G., Bopp, L., Canadell, J. G., Chevallier, F., Ciais, P., Ellis, R., Gloor, M., Peylin, P., Piao, S. L., Le Quéré, C., Smith, B., Zhu, Z., and Myneni, R.: Recent trends and drivers of regional sources and sinks of carbon dioxide, Biogeosciences, 12, 653–679, https://doi.org/10.5194/bg-12-653-2015, 2015. a

Smith, B., Prentice, I. C., and Sykes, M.: Representation of vegetation dynamics in modelling terrestrial ecosystems: comparison two contrasting approaches within European climate space, Global Ecol. Biogeogr., 10, 621–637, https://doi.org/10.1046/j.1466-822X.2001.t01-1-00256.x, 2001. a, b

Smith, B., Wårlind, D., Arneth, A., Hickler, T., Leadley, P., Siltberg, J., and Zaehle, S.: Implications of incorporating N cycling and N limitations on primary production in an individual-based dynamic vegetation model, Biogeosciences, 11, 2027–2054, https://doi.org/10.5194/bg-11-2027-2014, 2014. a, b

Soo‐Hyung, K. and Lieth, J. H.: A Coupled Model of Photosynthesis, Stomatal Conductance and Transpiration for a Rose Leaf (*Rosa hybrida* L.), Annals of Botany, 91, 771–781, https://doi.org/10.1093/aob/mcg080, 2003. a

van der Werf, G. R., Randerson, J. T., Collatz, G. J., Giglio, L., Kasibhatla, P. S., Arellano, A. F., Olsen, S. C., and Kasischke, E. S.: Continental-Scale Partitioning of Fire Emissions During the 1997 to 2001 El Niño/La Niña Period, Science, 303, 73–76, https://doi.org/10.1126/science.1090753, 2004. a

Vitousek, P. M., Ehrlich, P. R., Ehrlich, A. H., and Matson, P. A.: Human Appropriation of the Products of Photosynthesis, BioScience, 36, 368–373, https://doi.org/10.2307/1310258, 1986. a

von Bloh , W., Rost, S., Gerten, D., and Lucht, W.: Efficient parallelization of a dynamic global vegetation model with river routing, Environ. Modell. Softw., 25, 685–690, https://doi.org/10.1016/j.envsoft.2009.11.012, 2010. a

Walker, A. P., Johnson, A. L., Rogers, A., Anderson, J., Bridges, R. A., Fisher, R. A., Lu, D., Ricciuto, D. M., Serbin, S. P., and Ye, M.: Multi-hypothesis comparison of Farquhar and Collatz photosynthesis models reveals the unexpected influence of empirical assumptions at leaf and global scales, Glob. Change Biol., 27, 804–822, https://doi.org/10.1111/gcb.15366, 2021. a, b

- Abstract
- Introduction
- Solution of nonlinear equations
- Application to the problem
- Numerical performance and discussion
- Conclusions
- Appendix A: Parameters in photosynthesis
- Appendix B: Programming
- Appendix C: LPJmL v4 benchmark results
- Code and data availability
- Author contributions
- Competing interests
- Disclaimer
- Acknowledgements
- Financial support
- Review statement
- References
- Supplement

- Abstract
- Introduction
- Solution of nonlinear equations
- Application to the problem
- Numerical performance and discussion
- Conclusions
- Appendix A: Parameters in photosynthesis
- Appendix B: Programming
- Appendix C: LPJmL v4 benchmark results
- Code and data availability
- Author contributions
- Competing interests
- Disclaimer
- Acknowledgements
- Financial support
- Review statement
- References
- Supplement