This manuscript describes a new stomatal model that does not assume steady state responses to climate. Instead, this approach models gs with prognostic updates at each time step, which allows the model simulate limitations to gas exchange caused by stomatal speed and simplify leaf energy balance calculations.

The model was coupled with an LSM and evaluated at leaf and ecosystem levels. At leaf level the model could predict better the stomatal responses to changes in light regime. At ecosystem level, the model makes only a small difference at the canopy fluxes at mornings and evenings.

The paper is well written, and the results are generally clearly presented. I find the results of the paper interesting, and I really like the idea of the dynamic stomatal model because it seems a better representation of plant stomatal behaviour in a (rapidly) changing environment.  I believe this approach have the potential to improve how LSMs represent ecosystem carbon and water dynamics. My only criticisms are relatively minor methodological issues explained in detail below. Additionally, I think the paper would be improved if the authors had also validated their dynamic model against ecosystem level observations (i.e. eddy flux data) similarly to what they did for leaf-level data.

Author response: We would like to express our gratitude to the reviewer for taking the time to review our manuscript and provide valuable feedback. We agree with the reviewer that further validation on site-level measurements would be an improvement to illustrate the potential of our dynamic model. It can be the focus of follow-up research when detailed site-level datasets of eddy fluxes, plant traits, and meteorological conditions are available for accurate parameter calibration. We will add this to the discussion section for future research directions.

Specific comments:

L9: In the results you claim the daily effect of the dynamic model is notable (L244)

Author response: Thanks for pointing this out. The overall effect is not significant, as the daily mean differences are mostly less than 2%, but depending on the variations of environment and when considering the mornings and afternoons, the impacts can be notable (eg. up to -7.4%). We acknowledge our current writing in these two lines can be confusing and will update the phrasing as follows for clarity.

“The differences in fluxes between the NSS and SS predictions were not significant when integrated over monthly periods … but are notable in certain diurnal cycles depending on the radiation and other conditions, especially when considering the sub-diurnal scale.”

L24: The idea of using optimization to predict stomatal behaviour goes back to Cowan & Farquhar 1977, so I am not sure if using the term “more recently” is really appropriated.

Author response: Thank you for the suggestion. The “more recently” in our writing mainly intended to refer to the implementation and usage of stomatal models in large-scale LSMs, as most LSMs have been using empirical models, such as CLM, LPJ-GUESS, etc. We agree when considering the overall history of the optimization approach, the term we used here is not the most appropriate one, we will replace this term with “Efforts have also been made to…” for clarity.

L25: Most of the models discussed in these papers do not really optimize water use efficiency; they optimize the balance between carbon gain with some penalty, which can be related to plant hydraulics, non-stomatal limitations to photosynthesis, etc.

Author response: Thank you for the correction. We acknowledge our current phrase is not a very accurate summary of the optimization approach. Optimal models optimize the trade-offs between carbon gain and a variety of penalties related to stomatal opening, which may not (only) include absolute water loss that WUE is calculated from. We will correct this phrase as follows for accuracy: “Efforts have also been made to constrain stomatal behavior from the principle of optimizing the trade-offs between carbon gain with the related penalty of stomatal opening.”

L34-35: I assume it depends when the “next change occurs”? It would be useful to provide some quantitative examples of the time scales of stomatal response to environmental changes to clarify that.

Author response: Thank you for the suggestion. As shown in the field radiation measurements (Figure 7 & 9), the light received by the canopy keeps changing from sunrise to sunset, especially with fluctuations during cloudy days. Meanwhile, the time constant of stomatal responses, based on both our results (Figure 3) and previous studies, can range from a few minutes to more than half an hour. This indicates stomata in the natural environment are often not in steady states. We will add these quantitative examples with references to our introduction.

L100: What steady-state model do you use to calculate the target gs?

Author response: We used the Ball-Berry and Medlyn model (L112). On the leaf level, we used the Ball-Berry model for the test of our dynamic model; at the canopy scale, we used the Medlyn model because of the availability of vegetation traits datasets in our study region. We realize our current writing may have caused confusion, and we will move this description from Section 2.1.2 to Section 2.1.1 for clarity. 

L120: How many leaves/plants were used?

Author response: We tested our model on a couple of leaves, and for the manuscript, we only included two leaves with distinct time constants to illustrate the concept of temporal responses and show the differences in stomata behaviors across leaves, as well as what this may lead to further variations in parameter estimation when using traditional methods with the steady-state assumption.

L140: This vcmax value should be micromol instead?  L205: check vcmax units

Author response: Thanks for the correction. We will fix the typos in these lines.

L170: Wouldn’t a linear interpolation homogenize the environmental conditions over time and “mask” the real importance of the non-steady state model? I assume the dynamic model would only make a bigger difference in environments with rapidly changing environmental conditions. Would that explain why you have a relatively small impact of the dynamic model on the leaf and canopy simulations?

Author response: We concur with the reviewer's observation that linear interpolations may homogenize variations and this could be a reason for the relatively small impacts. On the other hand, light fluctuations (for which we used high-resolution measurements) tend to be the most rapidly-changing environmental condition compared to other variables (e.g. VPD, soil moisture, etc.) that we applied linear interpolations. Thus, the effect of our interpolation may not be relatively significant. The magnitude of such impact can be assessed with measurements of higher temporal resolution for other env variables in future studies. We will also update the related methods section to clarify this point.

L175: This section needs more details on the LSM configuration for these simulations to allow reproducibility. For example, how canopy light diffusion, carbon allocation/vegetation dynamics and soil moisture were handled. Maybe add this information as supplementary material.

Author response: Thanks for the suggestion. More detailed information on the LSM configuration can be found in the references of the vegetation module of CliMA Land: Wang et al. (2021, 2023). Briefly, in this study, canopy radiative transfer is modeled with a vertically layered canopy scheme adapted from the Soil Canopy Observation of Photosynthesis and Energy fluxes model (SCOPE), which includes leaf angular distribution to simulate hyperspectral reflectance and transmittance. In the current version of CliMA Land we employed, vegetation properties are prescribed with global datasets, such as LAI and Vcmax (for ‘carbon allocation/vegetation dynamics’). Soil moisture is prescribed with ERA5 reanalysis data and a tuning factor is employed to link stomatal responses to soil moisture status. We will also summarize this information and include it in the supplementary material in the final manuscript.

L245: Do the daily differences “disappear” at monthly scales because you have opposite RD in different days that cancel each other out? If that is the case, it would still be interesting to show the total monthly differences between SS and NSS (using for example the mean absolute error between SS and NSS). I believe the differences between models throughout the month could still result in different trajectories for the vegetation over time, for example, if the higher carbon gain of a model in a given day resulted in more leaves being produced on that day, this model productivity advantage that would accumulate over the other one.

Author response: The reviewer is correct, in our case, some of the daily differences are opposite and the aggregation makes the monthly differences less significant. We also agree with the reviewer that temporal responses of stomata could result in accumulated effects on vegetation growth. As the version of CliMA Land we employed has not implemented vegetation dynamics from net carbon gain from daily fluxes, we focused on evaluating the impacts of stomatal response on the gas exchange in this study. We appreciate the reviewer’s suggestion and will extend the discussion section to include this as a direction for future research.

L276: I could not find the computational cost differences in the results.

Author response: At each time step in traditional simulations, iterations are needed for steady solutions. In CLM4, the default setting is a 3-step fix-point iteration, which would be at a similar computational cost to our dynamic model if run at 10min resolution (3 sub-steps within 30min). In Section 3.3, we demonstrated that the dynamic model can be stably run at a resolution of 10min (L268-270). Furthermore, studies (Sun et al., 2012) have suggested this default setting in CLM4 does not always converge, so additional cost or algorithm improvements are needed. This algorithm is updated in CLM4.5 and CLM5, but iterations are still required as the steady-state assumption remains. The nested loop (Figure 1a) can take up to 40 iterations at a single time step (Bonan, et al., 2018). Thus, we conclude that our dynamic modeling represents a comparable efficiency to traditional SS simulations while providing predictions at a finer resolution and eventually also facilitates prognostic leaf temperature treatment.

L295: It could be interesting to see the cumulative effect of this unnecessary water loss on xylem hydraulic damage in future studies.

Author response: Thank you for the suggestion. We also believe this could be an interesting and valuable direction for future research with LSMs that implemented parameterization of such effects. We will include this in the discussion of future studies.

Fig. 2: I don’t consider this figure really essential for the manuscript. Considering you already have 9 figures maybe it would be best to leave it as supplementary material.

Author response: Thanks for the suggestion, we agree with the reviewer that it would be better to have it in supplementary materials. We will move this figure to the supplementary file in the final manuscript.

Fig. 6: Its hard to visualize the RD, please use darker shading.

Author response: Thank you for the suggestion. We have updated the shading to a darker one for better visualization.

The manuscript by Liu et al. presents how including dynamic stomatal responses impacts gas exchange predictions from leaf to canopy. The topic is of high importance as research carried out at the leaf level has shown the importance of stomatal dynamics but there is still a debate about how these dynamics influence gas exchange at a larger scale. The manuscript is well written and I like that the authors are addressing this gap in the literature, however, I have some concerns about the methods used.

Author response: We appreciate the time and effort the reviewer has dedicated to reviewing our manuscript and the valuable comments they provided for improving our study. We understand the reviewer’s concern regarding our methods, and we appreciate the opportunity to explain and discuss these issues. We acknowledge that some concerns may have arisen due to the level of detail in our methods section. Thus, we would like to provide more detailed descriptions and explanations of the methods we used in this study here, and we will also revise our methods section and add supplementary materials on our model configuration for clarity in the final manuscript.

What should be clearer is what are the challenges to simulating gas exchange at different scales from leaf to canopy and from seconds to months. For example, during a diurnal period, the slow temporal response of stomatal conductance will lead to the transient limitation of photosynthesis and potential damage to the PSII reaction centres (due to excessive energy received compared to the sink strength). These effects accumulate during the day and impact plant growth, which will in turn affect the gas exchange of the following days. These effects are not clearly accounted for by the model. Averaging the fluxes over a long period is useful for quantifying the ecosystem exchange but does not reflect that every day the environmental conditions will impact the leaf functioning and acclimation.

Author response: We acknowledge the valid point raised by the reviewer regarding the limitations of our model. The slow temporal response of stomatal conductance can have impacts on plants in various aspects, which, as the reviewer mentioned, include the direct gas exchange, as well as the potential damage to the PSII reaction centers and other accumulated effects on leaf functioning. Here, our focus is on taking the first step to scale up the stomatal response effects on the gas exchange from leaf to canopy and ecosystem in an LSM. We agree that when accurate parameterization of other effects is available, future efforts can be made to quantify those impacts on the canopy scale and provide further understanding of how dynamic stomata response affects vegetation dynamics in the long term. We will extend our discussion section and address these challenges and limitations, as well as suggest future efforts in these directions.

The effect of stomatal dynamics is also dependent on the fluctuation of the environment. The more it fluctuates, the more the stomatal behaviour will matter in the gas exchange.

Author response: We agree with the reviewer that the relative effects of stomatal dynamics depend on the fluctuation of the environment, as this is the main reason we employed high-resolution measurements of radiation for the comparison simulations. We will improve our clarity in addressing this point in the methods section.

The methods section is missing important information and requires more clarity.

Author response: We acknowledge the concern that the reviewer brought out about the level of detail in our methods section, and we would like to clarify unclear descriptions the reviewer has addressed point by point (as in the following responses). We will also summarize and add these to the methods section or supplementary materials in the revised manuscript.

From lines 90 to 105, it is not clear how gbc was set or calculated.

Author response: For leaf-level predictions, gbc is prescribed with the estimated gbc from LI6800 measurements. For canopy-scale simulations, in the CliMA Land version we employed in this study, we used a reasonable and fixed gbc value of 3/1.35 (gbw is assumed to be 3), which is relatively high conductance to make sure the gbc is not the main limiting factor of CO₂ supply (as our focus is on effects from stomatal conductance). We acknowledge more realistic gbc including calculation from wind speed and leaf width, which can be challenging on the canopy scale. As the vegetation module of CliMA Land implemented a vertically layered canopy scheme (rather than a big-leaf scheme), such calculation would require resolved wind speed at each layer and coupling with the atmospheric module, which is still under development. Thus, we set a reasonable gbc for all layers, which we believe will not result in significant impacts on our results, as the focus is the comparison of the two modeling approaches, rather than the absolute fluxes.

How the simplified model was solved and with which procedure is missing. Did the authors use an ODE solver?

Author response: The simplified model is solved with the Euler method with a fixed step size (the time steps used in simulations), as suggested in Figure 1. We will add a description to the supplementary material in the final manuscript, and also provide an example script of how we run the simplified model on the leaf level for clarity, along with other scripts for generating the results in this study.

How are the environmental variables included in this model? Did the authors use them as forcing variables that change continuously?

Author response: Yes, environmental variables from ERA5 reanalysis dataset (e.g. air temperature, dew-point temperature, volumetric soil water, wind speed etc.) are used as meteorological drivers and updated accordingly at each time step. We will clarify in the revised manuscript.

How is the steady state gs calculated? Which model was used: Ball or Medlyn? How were they calibrated?

Author response: We used the Ball-Berry and Medlyn model (L112). On the leaf level, we used the Ball-Berry model for the test of our dynamic model, and calibrated with the framework we described in Section 2.2.2 (L125-146). At the canopy scale, we used the Medlyn model because of the availability of vegetation traits datasets in our study region. We realize our current writing may have caused confusion, and we will move this description from Section 2.1.2 to Section 2.1.1.

At different places in the manuscript, the time steps used are different and without a clear explanation of why, it is confusing.

Author response: Thanks for pointing this out. In our leaf-level runs: as the length of total measurements is around 4.5 hours, to illustrate the effects of dynamic stomatal conductance response with higher accuracy, we used a relatively fine resolution of 10s.

On the canopy scale, for the comparison of differences in gas exchange, considering the balance between computational cost and accuracy, we chose a practical resolution of 1min for non-steady-state runs and compared it with traditional steady-state simulations at 30min resolution (the commonly used time step of LSMs). For the comparison of model efficiency, we also ran our model on 2, 6, and 10min (Figure 11) resolution to test model stability.

The authors did not include leaf energy balance to calculate leaf temperature, but this can have an important impact on their simulations. The vapour pressure gradient between the leaf and atmosphere depends on this and drives the water exchange. This should be acknowledged at least as a limitation in the discussion. Overall, it is difficult to judge the validity of the simulations, because important details are missing.

Author response: Thank you for the suggestion. We agree with the reviewer the calculation of leaf temperature from energy balance impacts our simulation, and the stomatal response contributes to the latent heat fluxes, as indicated in Figure 1 and the discussion section (L319-325). The implementation of leaf energy balance for leaf temperature with dynamic stomatal response would also require a non-steady-state energy balance modeling (in traditional LSMs, this is done with nested loops for steady solutions, as in Figure 1), which is under development in the latest version CliMA Land. We also plan to include dynamic leaf temperature in future simulations, so that the leaf flux calculation in our LSM can become fully prognostic, and we will also be able to better quantify the effects of stomatal response. We will further clarify this point in the related part of our discussion section.

We hope our explanations above can provide the details that the reviewer has suggested should be included. We will also summarize the key LSM configuration info and add it as supplementary material in the final manuscript.

For the exponential model used here, there is a time constant for stomatal opening and another one for stomatal closure. Stomata have very different rapidity of response for opening and closing. Here the authors only used one time constant and only estimated it for the opening part. The asymmetry of the time constants is strongly influencing the hysteresis between the morning and afternoon part of the simulation. This should be corrected.

Author response: Thank you for the suggestion. We acknowledge the differences in response speed of stomata opening and closing (L314-315) and agree this may influence the hysteresis pattern between the morning and afternoon. However, (a) although applying different time constants can be easily done in our model, it can be hard to assume valid time constants for opening and closure on a canopy scale without sufficient measurements on the site level for parameter estimate (as the magnitude of hysteresis will be dependent on the relative differences between the two time constants). (b) the focus of our study is mainly to take the first step of implementing a dynamic stomatal model in LSMs and illustrate its implications. Our point in this part (as one of the implications) is to suggest the temporal response of stomata could be one of the explanations for the observed hysteresis pattern, as previous studies have often merely focused on the meteorological variables to explain this phenomenon and neglected the stomatal response effects. We believe further and thorough studies are needed to investigate whether or how much the temporal response of stomata can contribute to explaining observed hysteresis, but we consider this is out of the scope of our current study.

We will address the influence of time constant on the hysteresis and suggest future efforts to further investigate the asymmetry of the time constants on the canopy scale, and how this may contribute to the observed hysteresis patterns of canopy fluxes.

The authors concluded that their optimisation procedure improves the parameter estimation compared to traditional linear fitting methods. It is not clearly explained in the text how these methods differ and what are the improvements. The need to reach steady state gs to estimate the parameter values is known. Using a coupled dynamic gs model and steady-state target model to fit the data of a light response curve has been done previously. In general, it is an interesting result but not the focus of this study.

Author response: Our conclusion about this procedure is mainly that it can provide a valuable alternative approach for parameter estimation. We successfully employed this approach to calibrate the parameters for our dynamic model and reproduced the leaf response curves. Compared to traditional linear fitting methods, (a) the Bayesian nonlinear inversion framework can optimize multiple parameters based on a joint fit of An and gs response curves (L125-126). In the traditional linear fits, Vcmax is often estimated with A/Ci response curve, and stomatal parameters are often estimated separately with light response curves. This approach can also be employed to estimate parameters with both A/Ci and light response curves. (b) As the reviewer mentioned, the linear regression method requires reaching an equilibrium at each environmental condition, which can be time-consuming. The bias of estimation with too short of a time step has also been discussed in previous studies (L284-286). In the meantime, our fitting with the dynamic model does not require steady states in principle, which can help reduce the time investment in parameter calibration. Our main goal of this study is to implement the dynamic model in our LSM and illustrate some of its implications. As the reviewer suggested, this part is not our main focus of this study but one of the implications we would like to point out, that implementation of the dynamic model can provide an alternative way of parameter estimate without the steady-state assumption.

Line 25: the optimisation theory does not optimize water use efficiency.

Author response: Thank you for the correction. The other reviewer has also pointed out our inaccurate description of optimization models. We acknowledge they mainly optimize the balance between carbon gain and a variety of potential penalty functions related to stomatal opening, which may not (only) include absolute water loss that WUE accounts for. We will correct this phrase as follows for accuracy: “Efforts have also been made to constrain stomatal behavior from the principle of optimizing the trade-offs between carbon gain with the related penalty of stomatal opening.”