Author’s Response to the Reviewer # 1 Comments

We appreciate the time and efforts by the editor and Reviewers in reviewing this manuscript and the valuable suggestions offered. We have attempted to address all issues raised by the Reviewers and hope that the revised manuscript can satisfy the Reviewer’s comments and journal’s requirements. The bold text indicates the comments proposed by Reviewers and the black text refers to our response to the comments.

In this paper, Yao et al. developed and implemented a new plant hydraulic architecture module “NHA” into ORCHIDEE-CAN based on soil-root-stem-leaf water transport continuum and the relationship between PLC and tree mortality. They compared the model performance of NHA model with two previous versions of the model to prove the efficacy of the new model in capturing the change of sap flow, soil moisture content, and GPP under drought events. They also evaluated model results against field measurements of leaf water potentials, biomass and mortality rates from a tropical lowland rainforest in eastern Amazonia. Their results show great potential of the NHA model to capture the drought-related tree biomass loss and mortality for tropical forests. The new model represents the state-of-the-art development of plant hydraulic model and will be of interest to the research community and readers of GMD. This paper is well written, and the results are nicely presented. I have some general comments as below.

For improvement, first, they should fit their new model into a broader field of mechanistic plant hydraulic models. They mentioned some previous work such as SPA model and Xu et al. (2016) but it’s still not very clear how they were motivated, how the new model was built on, and what are the strengths and weaknesses of their new model compared with other similar plant hydraulic models. They had some discussion starting from Line 547, but adding more details would be great.

[Response] Thanks for your suggestions. Our new model mainly focus on the extension of plant hydraulics to the hydraulic failure induced tree mortality. We add comparison with other plant hydraulic models in format of table as Table A1. Our model solves a complete water potential profile including nodes of soil, soil in root zone, root, stem and leaf as well as the hydraulic conductance dynamics. One of our strong strength is the extension of plant hydraulic structure to the modeling of drought-induced tree mortality using continuous higher loss of stem conductance as the indicator of occurrence of tree mortality event. Admittedly, weakness does exist in our model, for example, parameter retrieval can be further realized through data assimilation that uses more benchmark. More optimization paradigm can be integrated into our model, which would benefit the parameterization process. Please see line 613-618 in the revised main text. 

Table A1 Plant hydraulics in major vegetation models. The column of validation indicates how the model performance be validated against observation.

Second, one of the key limitations of the usage of such plant hydraulic models is numerous parameters, as shown in Table 1 in this paper. The authors focused on one site simulation with well-recorded plants’ traits. However some topics such as how sensitive and uncertain these parameters are, and how to parameterize the model at the regional and global scales might be interesting to add to the discussion. The authors may find this paper relevant to their discussion:

Liu, Y., Kumar, M., Katul, G. G., Feng, X., & Konings, A. G. (2020). Plant hydraulics accentuates the effect of atmospheric moisture stress on transpiration. Nature Climate Change, 10(7), 691-695.

[Response] We acknowledge that the parameters corresponding to one site simulation are subject to large uncertainty to some extent. Parameters like Ψ₅₀ reflects the vulnerability of tree species, with less negative value meaning higher vulnerability to water stress. We made sensitivity test by varying degree of vulnerability (e.g. Ψgs₅₀). Figure R1 shows that when Ψgs₅₀ equals -1.2 MPa, the annual mortality rate would be more comparable with the observation.

We also add discussion regarding the parameterization of model at the regional and global scales. Generally, three means can be resorted to benefit such realizations. The first one can be embedding the plant trait database like TRY (Kattge et al., 2020) into our process-based model although the records are still limited in aspect of hydraulic traits. The second solution can be the optimization of hydraulic parameters using e.g. Monte Carlo Markov Chain with measurements or remote sensing products as constraints like traits retrieval in Liu et al. (2021b) or other data-assimilation system like ORCHIDAS. The third method can be building simple regression formula between plant traits and the climatology where the plants reside in. In next step, these solutions will be attempted to test the generalization of process-based model performance at large scale. Please see line 683-690 in the revised main text.

Figure R1 Modeled stem mortality rate with regard to different Ψgs_,50 values. when Ψgs₅₀ equals -1.2 MPa, the annual mortality rate would be more comparable with the observation.

Third, some references when the authors described the equations in Methods are missing.

[Response] We check the equations carefully and add the necessary references.

More information about throughfall exclusion experiment and model simulation set up is needed as well. Below, I provide more specific comments:

[Response] We add description for throughfall exclusion experiment and model simulation setup.

-Line 225: Any references for the sigmoidal relationship? How about other relationships such as linear, logistic, or exponential?

[Response] We add reference for the sigmoidal relationship (Pammenter and Van der Willigen, 1998). Linear may be not appropriate as we need impose segmentation on top of linear type to avoid the unrealistic values. Exponential and logistic functions can also be alternative since they can correspond to quick loss of conductance after reaching a critical threshold as well.

-Line 275: Please provide reference and a simple description for the gs model. L is not defined either.

[Response] The aim of this g_s model is to let g_s vary following dynamics of leaf water potential in sigmoidal function then g_s can be coupled into the plant water transport system via the transpiration supply. Meanwhile, the g_s is assured to close to 0 in the night, mediated by the radiation-related variable (.

g_s, g_max and g_min are in unit of mmol m^-2 s^-1.  is the function of short-wave radiation (Rad), which is used to ensure that g_s at night is very low. L and L_k are parameters specifying the strength of short-wave radiation limitation on stomatal conductance.

-Line 280: What’s the gs model in the SPA model, is that the same one used in this study?

[Response] The g_s model in the SPA model is different from what we used in this study. In SPA model, g_s in each canopy layer is obtained by maximizing the marginal carbon gain of stomatal openness, that is the optimization of intrinsic water use efficiency. Specifically, g_s in SPA increases iteratively until there is only negligible increase in assimilation or the leaf water potential is below a minimum threshold (Williams et al., 1996). We briefly mention the gs in SPA model in line 299-300 in the revised main text.

-Line 332: How is LAI modeled in this study?

[Response] LAI is determined by leaf mass, which is regulated by the leaf growth, leaf turnover and leaf loss due to drought-induced tree mortality. Please see new line 355-356 in the revised main text.

-Line 346: More information such as the plot size and duration of the experiment about the TFE site could be added here so readers don’t need to read the cited papers.

[Response] We add more description for the experiment. There are two experiments, which were carried out since the beginning of 2001. A throughfall exclusion experiment (TFE) started in the end of dry season in 2001, where 50% of canopy throughfall is excluded by plastic roof at the height of 1-2m above the ground. It is of 1-ha size. Another 1-ha control plot is also set without any manipulation. Here the observation data we used extends to 2008 at most due to data access issue, but these experiments are still running. 

-Line 353: What are the similarities and differences between SPA model and your model?

[Response] The similarities lie in the framework of hydrodynamic simulation, in which both models follow supply-demand theory. Specifically, the transpiration is limited by water supply derived from plant hydraulics. FvCB model is used for photosynthesis calculation in both models. The difference mainly relates to the computation of stomatal conductance (g_s) and water potential. SPA model represents g_s using optimization theory, while ORCHIDEE-CAN-NHA lets g_s vary with leaf water potential. The SPA model did not realize more explicit hydraulic segmentation, like root water potential was not modeled yet. As ORCHIDEE is a complex land surface model, including carbon allocation, phenology, turnover, tree mortality, SPA model did not resolve these processes yet.

-Line 360: What meteorological forcing was used to drive the model, at what temporal resolution? Were the simulations coupled with a climate model or offline?

[Response] The meteorological forcing is of half-hourly time step. The half-hourly meteorological data are measured using an automatic weather station located at the top (51.5 m) of a tower 1 km from the experimental plot. The simulation was ran offline without coupling with a climate model.

How was the TFE simulation carried out? Was the precipitation be reduced to 50% of CTL level at each model time step?

[Response] Yes. In model simulation, for TFE setup, we modify the precipitation forcing by cutting 50% of precipitation at each half-hourly time step while keeping temperature and down-ward short-wave radiation unchanged.

Is the model also initialized with real forest inventory data? How do 20 circumference classes correspond to the real-world situation?

[Response] The model was not initialized with real forest inventory data. In model routine, the tree density in each size class is firstly prescribed as a function related to PFT-related maximum tree height and initial number of trees. Then the tree density in each class changes with carbon allocation during tree growth, i.e., trees would move from current class to the next one. These 20 circumference classes constitute the demographic structure, spanning from smaller trees to bigger one. As the existing studies at Caxiuana site did not report its demographic structure, we retrieve the observed tree density distribution from annual mortality rate in each size group and total annual mortality rate. We assume here F (fraction of each size group) did not show inter-annual variation. This is an approximate estimation, from which the sum of F1 to F3 could not equal 1 due to the above assumption.

F₁m₁₁+F₂m₂₁+F₃m₃₁=T₁

F₁m₁₂+F₂m₂₂+F₃m₃₂=T₂

F₁m₁₃+F₂m₂₃+F₃m₃₃=T₃

F_j is the fraction of group j (j=1,2,3). T_i is annual tree mortality rate in year i. Three groups correspond to class with DBH < 20cm, class with DBH between 20cm and 40cm, and class with DBH above 40cm.

a_ji is the annual mortality rate of group j in year i (i=1,2,3). Here we choose mortality rate in three years to solve F₁, F₂, and F₃ as only square matrix (here 3 rows and 3 columns) has the inverse matrix.

F could vary according to different T as we can choose three-year data. Generally, F₁ is the highest and F₃ is the lowest. For example, when we choose year 2002 to 2004, the fraction of first group is 62% (F₁ = 0.62, F₂ = 0.25, F₃ = 0.07). In model output, the group with DBH less than 20cm accounts for almost 70% of all tree individuals, which is similar to the real situation at Caxiuana site. Modeled fraction of each size group also decreases with tree size (Figure R2).

Figure R2 Tree density distribution of 20 circumference classes.

-Line 365: Past tense for “run”, and “compare” in line 367.

[Response] They have been revised. Thanks.

-Line 384: Could the authors discuss why their new model underestimated sap flow in the dry season but overestimated it under TFE conditions?

[Response] The simulated transpiration could be limited by water supply (water limitation) or water demand (energy limitation). Under CTL, there is almost no water limitation even in dry season. The underestimated sap flow can be due to that model tends to underestimate the sensitivity to VPD increase in dry season. Under TFE, there is water supply limitation. The possible reasons for such overestimation under TFE can be that the sensitivity of water supply to drop in soil moisture is underestimated or the too slow soil water drainage in our model setup relative to that in reality (Kennedy et al., 2019). Please see line 420-425 in the revised main text.  

-Line 421: What mechanism leads to the larger seasonal amplitude of modeled GPP compared with SPA model?

[Response] In SPA model, GPP is simulated using FvCB model regulated by optimization of intrinsic water use efficiency, in which the optimization target is  (A is assimilation, g_s is stomatal conductance), not accounting for VPD. So the magnitude of GPP variation would not be too high. In ORCHIDEE-CAN-NHA that we used here, larger seasonal amplitude of modeled GPP especially the low GPP in dry season under TFE is due to higher water limitation arised from our hydraulic architecture. Please see the explanation in line 475-479 in the revised main text.

-Line 550: What’s the leaf-level demand of Xu et al. (2016)?

[Response] Leaf-level demand of Xu et al (2016) is transpiration. It is calculated as the product of total conductance for water vapor, total leaf area per cohort and the gradient of water vapor concentration between leaf intercellular space and leaf boundary layer.

-Figure 2: Color for Ψ50 = -1.6 is too weak to be seen.

[Response] We revise the Figure 2 to let each line be clear now. We also move Figure 2 to SI following reviewer 2’s comment.  

Figure R3 Sigmoidal relationship between stem sapwood conductance (normalized by total leaf area) and stem water potential. Line colors correspond to different Ψ₅₀ values. The line types (continuous, dashed and dotted line) denote different curvature parameters (a_stem).

References

Christoffersen, B. O., Gloor, M., Fauset, S., Fyllas, N. M., Galbraith, D. R., Baker, T. R., Kruijt, B., Rowland, L., Fisher, R. A., and Binks, O. J.: Linking hydraulic traits to tropical forest function in a size-structured and trait-driven model (TFS v. 1-Hydro), Geoscientific Model Development, 9, 4227-4255, 2016.

Cochard, H., Pimont, F., Ruffault, J., and Martin-StPaul, N.: SurEau: a mechanistic model of plant water relations under extreme drought, Annals of Forest Science, 78, 1-23, 2021.

De Kauwe, M. G., Medlyn, B. E., Ukkola, A. M., Mu, M., Sabot, M. E., Pitman, A. J., Meir, P., Cernusak, L., Rifai, S. W., and Choat, B.: Identifying areas at risk of drought‐induced tree mortality across South‐Eastern Australia, Global Change Biology, 2020.

Eller, C. B., Rowland, L., Mencuccini, M., Rosas, T., Williams, K., Harper, A., Medlyn, B. E., Wagner, Y., Klein, T., and Teodoro, G. S.: Stomatal optimization based on xylem hydraulics (SOX) improves land surface model simulation of vegetation responses to climate, New Phytologist, 226, 1622-1637, 2020.

Fisher, R., Williams, M., Da Costa, A. L., Malhi, Y., Da Costa, R., Almeida, S., and Meir, P.: The response of an Eastern Amazonian rain forest to drought stress: results and modelling analyses from a throughfall exclusion experiment, Global Change Biology, 13, 2361-2378, 2007.

Kattge, J., Bönisch, G., Díaz, S., Lavorel, S., Prentice, I. C., Leadley, P., Tautenhahn, S., Werner, G. D., Aakala, T., and Abedi, M.: TRY plant trait database–enhanced coverage and open access, Global change biology, 26, 119-188, 2020.

Kennedy, D., Swenson, S., Oleson, K. W., Lawrence, D. M., Fisher, R., Lola da Costa, A. C., and Gentine, P.: Implementing plant hydraulics in the community land model, version 5, Journal of Advances in Modeling Earth Systems, 11, 485-513, 2019.

Liu, Q., Peng, C., Schneider, R., Cyr, D., Liu, Z., Zhou, X., and Kneeshaw, D.: TRIPLEX-Mortality model for simulating drought-induced tree mortality in boreal forests: Model development and evaluation, Ecological Modelling, 455, 109652, 2021a.

Liu, Y., Parolari, A. J., Kumar, M., Huang, C.-W., Katul, G. G., and Porporato, A.: Increasing atmospheric humidity and CO2 concentration alleviate forest mortality risk, Proceedings of the National Academy of Sciences, 114, 9918-9923, 2017.

Liu, Y., Holtzman, N. M., and Konings, A. G.: Global ecosystem-scale plant hydraulic traits retrieved using model–data fusion, Hydrology and Earth System Sciences, 25, 2399-2417, 2021b.

Naudts, K., Ryder, J., McGrath, M., Otto, J., Chen, Y., Valade, A., Bellasen, V., Berhongaray, G., Bönisch, G., and Campioli, M.: A vertically discretised canopy description for ORCHIDEE (SVN r2290) and the modifications to the energy, water and carbon fluxes, Geoscientific Model Development, 8, 2035-2065, 2015.

Pammenter, N. v., and Van der Willigen, C.: A mathematical and statistical analysis of the curves illustrating vulnerability of xylem to cavitation, Tree physiology, 18, 589-593, 1998.

Papastefanou, P., Zang, C. S., Pugh, T. A., Liu, D., Grams, T. E., Hickler, T., and Rammig, A.: A dynamic model for strategies and dynamics of plant water-potential regulation under drought conditions, Frontiers in Plant Science, 11, 373, 2020.

Wang, Y., Köhler, P., He, L., Doughty, R., Braghiere, R. K., Wood, J. D., and Frankenberg, C.: Testing stomatal models at stand level in deciduous angiosperm and evergreen gymnosperm forests using CliMA Land (v0. 1), Geoscientific Model Development Discussions, 1-35, 2021.

Williams, M., Rastetter, E., Fernandes, D., Goulden, M., Wofsy, S., Shaver, G., Melillo, J., Munger, J., Fan, S. M., and Nadelhoffer, K.: Modelling the soil‐plant‐atmosphere continuum in a Quercus–Acer stand at Harvard Forest: the regulation of stomatal conductance by light, nitrogen and soil/plant hydraulic properties, Plant, Cell & Environment, 19, 911-927, 1996.

Xu, X., Medvigy, D., Powers, J. S., Becknell, J. M., and Guan, K.: Diversity in plant hydraulic traits explains seasonal and inter‐annual variations of vegetation dynamics in seasonally dry tropical forests, New Phytologist, 212, 80-95, 2016.

­­

Author’s Response to the Reviewer #2 Comments

We appreciate the time and efforts by the editor and Reviewers in reviewing this manuscript and the valuable suggestions offered. We have attempted to address all issues raised by the Reviewers and hope that the revised manuscript can satisfy the Reviewer’s comments and journal’s requirements. The bold text indicates the comments proposed by Reviewers and the regular text refers to our response to the comments.

This paper describes the implementation of a hydraulics scheme into the land surface model ORCHIDEE and its evaluation against the Caxiuana drought experiment. In general, this is a nice piece of work, but the presentation could be quite significantly improved.

One important point is the need to show how the model simulations have changed since the new routines were added. The figures only show output from the new version of the model. To assess the value of the added subroutines, the paper needs to show output from previous versions of the model for comparison. There is some in the supplementary, but there is insufficient quantitative assessment of how each version of the model performs. The R values for sap flow are lower in the new model version than the previous one, suggesting a degradation of model performance. The comparison of GPP with previous models is qualitative only. It would be valuable to add some statistics to compare performance of different model versions.

[Response] Thanks for the comments. We add comparison of transpiration as well as biomass loss with ORCHIDEE-CAN and ORCHIDEE-CAN-RS. We do not put water potential here as the other two model versions do not have such outputs. Here we show the updated results after increasing interception reservoir (following later comment from the Reviewer). We do not aim to let the model outputs perfectly match the observation in order to avoid the overfit issue. We add below evaluation of model performance in supplementary.

In terms of comparison on transpiration (Table R1), under CTL, the correlation coefficient with the observation is similar among three model versions (0.71-0.76) although there is indeed a bit increase in other error metrics in ORCHIDEE-CAN-NHA like root mean square error (RMSE) and mean absolute percentage error (MAPE). ORCHIDEE-CAN-NHA performs better in water stress condition (under TFE) in aspects of error metrics but shows a bit lower correlation with observation than other two versions.

As SPA-GPP is also based on model simulation, here we compare the seasonal magnitude of GPP among different model versions more specifically (Table R2), without using correlation and other errors metrics. During 2001-2003, under CTL, SPA-GPP shows very little difference between wet and dry season, while GPP from ORCHIDEE-CAN-NHA presents 1.1 gC m^-2 day^-1 difference, which is similar with two previous versions. Under TFE, SPA-GPP drops from wet to dry season, while ORCHIDEE-CAN-NHA shows small GPP increase in dry season as water stress effects tend to occur later in dry season, rather than the beginning of dry season like SPA. Dry season GPP increase is also found in other two model versions in spite of a bit difference in the magnitude.  

With regard to the biomass loss, only ORCHIDEE-CAN-NHA can produce a comparable amount of biomass loss under TFE by the end of 2008 (Table R3). 

Table R1 Model performance of transpiration evaluated by R, RMSE (unit: mm d^-1), and MAPE against observed transpiration.

Table R2 Comparison of modelled GPP from 2001 to 2003. Since SPA-GPP is based on model simulation as well, here we just compare the magnitude of GPP, without using correlation and other errors metrics. Unit : gC m^-2 day^-1.

Table R3 Model performance in simulation of biomass loss. Here the values in the table refer to relative change to the beginning of the experiment (2001).

There were quite a few questions about the model description. It would be of great value to go through the symbols used and try to make them consistent, instead of using a mixture of abbreviations and symbols. It is confusing to have WD the wood density, rho-root the root density, WC the amount of water per unit volume of sapwood, and rootwc the amount of water per gram root biomass. Try to come up with a more systematic set of symbols. In particular, avoid abbreviations instead of symbols (e.g. use D rather than dbh) and avoid using variable names from code such as circ_class_mor or counterPLC50. Give these symbols. Also use capitals consistently, e.g. Cleaf, Cstem and Croot should all have capital “C”

[Response] Thanks for the comments. We revised the model description part by using symbols especially those variables whose previous name is its abbreviation.

Specific leaf area: S, Leaf dry matter content: L, mass of water per unit of sapwood volume: γ, wood density: δ, root shoot ratio: θ, root water content: ε

Ensure to give all units clearly in text and ensure they are consistent. For example, is capacitance in units of mmol (line 203) or mmol m-2 MPa-1 (line 207) or in kg m-3 MPa-1 (Table A1)? I suggest checking over all of the equations thoroughly to ensure units are correct throughout the text.

[Response] Thanks for the reminder. The unit mentioned in line 203 is the unit of water storage. To avoid ambiguity, I specify the unit more clearly. The unit of stem and root capacitance is written as kg m^-3 MPa^-1. In the model, we also did a unit transform from kg to mmol. I specify it in the description of model formula.  

Eqn 4: m_sap,max = v_stem *WC

Vstem is the volume of a cylinder of diameter DBH and height h (eqn 6) so overestimates volume of a stem. How is the stem form factor corrected for? How is this then converted to sapwood?

[Response] We acknowledge that using cylinder as approximation for the tree volume will cause overestimation. The amount of water per unit stem volume is from Suzuki (1999), which is defined as . 

WC is defined as the mass of water per unit sapwood volume in mol m-3. It should be defined as the maximum mass of water, or the mass of water when water potential = 0. Clarify, is this per unit sapwood volume or per unit stem volume?

[Response] WC is from Suzuki (1999), which is defined as. We correct that WC (now symbol γ) corresponds to the maximum mass of water per stem volume.

Figure 2 does not seem important or relevant enough to include as a main figure. It just shows the form of the sigmoidal relationship for different parameter values. The different values are not used in the paper, however, so it’s not clear why this wide range of parameter values are shown.

[Response] We aim to give an illustration of the sigmoidal relationship between hydraulic conductance and water potential as well as how this relationship varies with shape parameter and Ψ₅₀. We moved Figure 2 to SI.

What happens when the canopy is wet? I note that in Figure 4, the canopy evaporation is a tiny fraction of ET, which seems very unlikely for this wet, high-LAI forest. These numbers need a reality check

[Response] When the canopy is wet, there are both leaf-level transpiration and canopy evaporation (from intercepted leaf water or dew re-evaporation – the model produces dew under certain conditions when the air is more humid than the surface from night to predawn). In our model, the potential fraction of rain that can be intercepted by the canopy is set to 30% for tropical evergreen forest. The interception is also regulated by the amount of interception reservoir, which is related to LAI and a coefficient λ (transforming LAI into size of interception reservoir, (Carlyle-Moses and Gash, 2011)). We find that in previous setup, the low canopy evaporation is mainly due to a very small interception reservoir (small λ). Measurement at Manaus site showed the interception loss is about 8.9% of annual rainfall (Lloyd et al., 1988). As we do not find other evidence at Caxiuana site, we use this Manaus value as reference. Here we increase the interception reservoir and now interception is about ~8% of annual rainfall (Figure R1). All the results are updated correspondingly.

Figure R1 Modeled (ORCHIDEE-CAN-NHA) daily soil evaporation (E), canopy evaporation (CE) and transpiration (T) during 2001-2008. The arrows point to the start of TFE in the beginning of 2002.   

I found the representation of Tdemand (eqn 23) to be remarkably simple – one would normally expect a land surface model such as Orchidee to have a more complex representation of T, including a boundary layer conductance and some scaling of gs to the canopy. Is the canopy transpiration the same as Tdemand?

[Response] Canopy transpiration is not equal to T_demand (potential transpiration demand, we rename this variable as PT_demand later for clarity) in our model. PT_demand acts as a constraint for the transpiration in water demand aspect. Regarding the transpiration, in ORCHIDEE, it is related to canopy surface resistance, vegetation structural resistance, air density, aerodynamic conductance and relative air humidity, thus state of the art parameterizations.   

Does this value of gs affect assimilation? How has the assimilation (and GPP) calculation changed?

[Response] This g_s value (varying with leaf water potential) does not affect assimilation directly but relates to the water supply for transpiration. Normally, when there is no water limitation, GPP is calculated by FvCB model that is the minimum of the Rubisco-limited rate of CO₂ assimilation and the electron transport-limited rate of CO₂ assimilation (Yin and Struik, 2009). If the water supply cannot meet the transpiration demand, then the transpiration would be re-calculated to match the water supply part. The canopy resistance is also re-calculated and so does the GPP.  

Be more specific about how water potentials are found.

[Response] Water potentials are solved to let the water supply equal to water demand at each organ. In the model, HYBRD1 function from Minpack package in Fortran is used, which seeks a zero of N nonlinear equations in N variables. The evaluated function is the difference between water supply and water demand at each organ level. This function iteratively minimizes the absolute value of the evaluated function. The initial estimate of the solution vector is quite important and comes from the water potential at last time step. For example, the initial estimate for leaf water potential at time step t that will be used in the formula is the stem water potential at time step t-1.

Line 270: we decrease leaf water potential until the difference between leaf water supply and demand is “close to zero” – How does this algorithm work? How close is tolerable?

[Response] Leaf water potential is solved using HYBRD1 function (see above response). The tolerance is 0.00001MPa. When the relative error between two consecutive iterates is below the tolerance, the calculation routine is terminated.

Line 294 and line 304: we “try to solve” Why only “try” ? How is the water potential found, and what happens if one can’t be found?

[Response] ‘try’ is removed from this sentence to avoid ambiguity. The water potential is solved using HYBRD1 routine that finds the solution for a series of non-linear equation. If no solution can meet current requirement and the number of iterations reaches the limit, water potential at final calculation step would be used. 

Give some indication of how parameter values are chosen. The table does list references, but it is not clear how values are chosen from the references.

[Response] We did some sensitivity tests by attempting different values combination of parameters within range of records in literatures, like degree of vulnerability, Ψ₅₀, and degree of sensitivity, a (shape parameter). Parameters set that can better capture the observed variation of drought-induced tree mortality (especially the higher tree mortality rate in larger cohorts) were chosen. Figure R2 shows that when Ψgs₅₀ equals -1.2 MPa, the annual mortality rate would be more comparable with the observation. We do not aim for a perfect match between model output and observation to avoid the overfit issue during the generalization of the model. 

Figure R2 Modeled tree mortality rate varies with different parameters setup (Ψgs₅₀ was taken as example here).

It’s unfortunately not acceptable to refer to other papers that are still in review. The Joetzjer et al. (in review) paper was not accepted in Biogeosciences in 2018.

The fact it has not yet appeared raises some questions. This paper does rely quite heavily on that one, so it seems essential that that paper be accepted before this one can be. There may of course be some extenuating circumstances.

[Response] Joetzjer et al (2022) has been published by Ecological Modelling at doi: 10.1016/j.ecolmodel.2022.109969. 

Section 2.1.6 is quite disconnected from the rest of the model implementation and it is not clear what has changed here from previous versions of the model.

[Response] This section describes the related carbon fluxes we used later in model-observation comparison. This calculation process is the same among different model versions and we put it in a separated section 2.2 now.

It would be valuable to add more interpretation of the outputs of the model in terms of underlying assumptions. For example, it’s noted that leaf water potentials are lower in the taller trees. The effect of height should be about -0.1MPa / 10m. Once this is accounted for there are similar LWPs across cohorts, which is somewhat surprising given that cohorts have different rooting depth and see different soil moisture. There also doesn’t appear to be a lot of difference in the PLC by cohort (Figure 9). The discussion later talks about the larger mortality rates in large trees, but it’s not clear how this arises from the model structure. It would be useful to talk through how this works in the text.

[Response] We assume that smaller cohorts can only access to shallower water table whereas larger cohorts can absorb soil water from the deeper soil layer. Soil water potential in root zone (Ψ_soil-root) is weighted by the amount of water that can be absorbed from each soil layer (ɳ). ɳ is determined by soil water potential and also the soil-root resistance. Ψ_soil-root did not show too much variation among different cohorts (Figure R3). Then the leaf water potential difference among cohorts is mainly contributed by the height effect. 

With regard to the tree mortality sub-model, we assume when the cumulated drought exposure reaches our specified threshold (15 continuous days with PLC above 50%), then there would be one mortality event. Such requirement is further relaxed by allowing 5 days wet break during a water stress condition, which would not impede the accumulation of the drought exposure. Figure R4 and R5 show that smaller cohort (#5 here) shows a bit larger variation in water potential dynamics and corresponding PLC, which indicates that an adequate cumulated drought exposure occurs less frequently than that of larger cohorts (#20 here). Thus the higher annual tree mortality rate is found in larger cohorts.

Figure R3 Soil water potential in root zone (Ψ_soil-root). The top panel shows the temporal distribution of Ψ_soil-root at daily time scale. The bottom one shows monthly Ψ_soil-root in cohort #5 and #20 for better comparison of their dynamics. 

Figure R4 Percentage of loss in stem hydraulic conductance (at daily time scale) and drought exposure in year 2005.  Cohort #5 and # 20 are shown.

Figure R5 The temporal dynamics of Ψ_leaf, Ψ_stem and Ψ_root in year 2005.  Cohort #5 and #20 are shown.

It does also seem odd that the lower soil layers dry out much more than the upper soil layers. It seems that the plants are preferentially using water from lower in the soil profile. Again, it would be useful to talk through what the model is doing in terms of water uptake.

[Response] The soil moisture content (SMC) at each layer is influenced by infiltration, evaporation, transpiration and drainage. The amount of water that can be absorbed from each layer (ɳ) is determined by its water potential and also soil-root resistance. Soil water potential decreases with soil depth while soil-root resistance becomes much smaller with soil depth as well. Therefore, ɳ does not change monotonically with soil depth. For example, in wet season in 2005 under TFE, ɳ in deeper soil layer is higher than that in top layer. While in dry season, ɳ in deeper soil layer can decrease to almost 0, when the water supply mainly comes from the shallower layer. In year 2004, even in dry season, lower soil layers can contribute a lot to water uptake (Figure R6). Soil water potential in root zone (Ψ_soil-root) is weighted sum of soil water potential in each layer and ɳ in each layer.

Figure R6 Distribution of ɳ in each soil layer (year 2004 and 2005 are taken as examples here). Top panels show ɳ in top three layers (0m, 0.002m 0.006m), and bottom panels show ɳ in lower layers (0.248m, 0.499m, 0.999m). The period with shading corresponds to the dry season from July to November. 

Smaller points:

Line 200 mentions “the first time-step” but is that just the very first half hour of a ten year simulation or is it every day? If water potentials are assumed the same in the first time step, what value do they take?

[Response] “the first time-step” points to the very first half-hour of the simulation. At first time-step, the initial value of Ψ_leaf, Ψ_stem, and Ψ_root all equal to Ψ_soil-root, which is the weighted sum of soil water potential. 

Please justify eqns 14 and 15.

   (14) 

  (15)

[Response] The value 2 in front of k_stem,t in each equation denotes that only half of stem is accounted for in upper part (water flow from middle of the stem to leaves) and trunk part separately (water flow from middle of the root to middle of the stem). Half of root length is considered in trunk part as well. The water transport process is assumed to be similar to electric current, of which the resistance (the reciprocal of hydraulic conductance) should be added up along the water transport path. 

Line 246: Please give correct units for J. (mmol m-2 s-1 ?)

[Response] We add the unit for J. The unit of J is mmol.

Please include values and units for the parameters in eqn 24.

[Response] The values and units for the parameters in eqn24 are shown in Table A2.

Line 349, “morality” should be “mortality”

[Response] Thanks.

Figure 3, how is sapflow extracted from the model? Is it the same as “T” in Figure 4 and “Tsupply” in Figure S4, or are these different outputs?

[Response] Here sap flow is treated as transpiration. The sap flow in Figure 3 is the same as ‘T’ in Figure 4. T_supply is transpiration supply, which is larger than transpiration when there is no water stress.

Please give full figure captions for the supplementary material. What is shown in Figure S1, exactly? What do the grey bands represent? Are the values given on a half-hourly or daily basis, and if daily, how are they averaged? Do the gs values differ by cohort?

[Response] We revise the figure captions in supplementary information. Here this g_s value differs by cohort as it is also regulated by the leaf water potential per cohort. I chose cohort #10 as an example here. The value is averaged to daily mean directly. The gray bands indicate dry season from July to November.  

What are the values shown in Figure S2? Were the data from Lin et al. filtered to show just high-PAR values? Are the observational values in fact comparable with the modelled values?

[Response] Figure S2 shows the g_s data from observation and our model output. The g_s data collected in Lin et al (2015) only provides the records during the time period we show here rather than the whole-year records. Lin et al (2015) did not mention the filtering work. Admittedly, there is model-observation mismatch in term of the magnitude of intra-annual variation but the magnitude of modelled values falls in the range reported from observation.

What are the measured and modelled values shown in Figure S8? How are the modelled values averaged over the cohorts? Where are the obs measured, and what is the uncertainty?

[Response] Figure S8 shows the leaf water potential values at Caxiuana site, which were monitored by a digital pressure bomb (Fisher et al., 2006). The measurement uncertainty can be due to limited sampling of leaves. The modelled value is an example from cohort#10. We did not average the value over all cohorts. 

I did go looking at the code repository but it’s very large and not clear where the new code resides. It would be useful to indicate which subroutines were modified / added in this version of the code

[Response] I add description mentioning the modification in this new model version. Please see the ‘details of model improvement’ section in below link:

https://forge.ipsl.jussieu.fr/orchidee/wiki/GroupActivities/CodeAvalaibilityPublication/ORCHIDEE_CAN_NHA

Details of model improvement:

The improvement of the new version mainly focuses on the simulation of dynamic water potential profile, hydraulic conductance, water storage dynamics regulated by capacitance, and drought-induced tree mortality in aspects of tree density and biomass.

Transpiration supply is updated in ./src_sechiba/sechiba_hydrol_arch.f90. Dynamic water potential profile and percentage loss of stem hydraulic conductance (PLC) are also solved in this script. Subroutine hydrol_arch is modified. Then PLC is transferred to./src_stomate/stomate_mark_kill.f90, where the number of continuous days with PLC above 50% is calculated as drought exposure. Then the amount of trees that would be killed is derived. The update of tree density and transfer of live biomass to litter are solved in ./src_stomate/stomate_kill.f90.   

References

Carlyle-Moses, D. E., and Gash, J. H.: Rainfall interception loss by forest canopies, in: Forest hydrology and biogeochemistry, Springer, 407-423, 2011.

Fisher, R. A., Williams, M., Do Vale, R. L., Da Costa, A. L., and Meir, P.: Evidence from Amazonian forests is consistent with isohydric control of leaf water potential, Plant, Cell & Environment, 29, 151-165, 2006.

Lloyd, C. R., Gash, J. H., and Shuttleworth, W. J.: The measurement and modelling of rainfall interception by Amazonian rain forest, Agricultural and Forest Meteorology, 43, 277-294, 1988.

Suzuki, E.: Diversity in specific gravity and water content of wood among Bornean tropical rainforest trees, Ecological Research, 14, 211-224, 1999.

Yin, X., and Struik, P.: C3 and C4 photosynthesis models: An overview from the perspective of crop modelling, NJAS-Wageningen Journal of Life Sciences, 57, 27-38, 2009.