General comments

The authors present a new hydrologic modeling platform over continental China, based on the ParFlow model, aimed at providing information for surface water and groundwater resources management. Their setup is adapted from the CONUS 2.0 modeling platform over the US. The authors discuss the parameters and input data used and they provide a comparison with modeled and observed data for groundwater table depth and river discharge.

The paper is well written and pleasant to read. The results are interesting and well presented via clear figures. My main concern is about the comparison with other datasets, which may be a little too simple, as detailed below.

In particular:

• l. 278-293: Do I understand correctly that for the spinup, the authors used the 1981-2010 P-ET average as constant atmospheric forcing until a quasi-steady state was reached? Is this resulting state used for the evaluation in the next section or have the authors simulated a transient run over 1981-2010, starting from this quasi-steady state? This is not clearly stated, while it is very important for the evaluation and its analysis. For the following comments, I will assume that the authors evaluate the resulting quasi-steady state against other data sets.
• section 4.1. and 4.2.: While the main motivations for this modeling platform are (1) the impacts on water resources of the increased frequency, intensity, and duration of extreme weather events and (2) the management of these water resources, e.g., to prevent water scarcity, the authors limit their evaluation to a comparison of the steady state, which represents an idealized situation that never happens in the real world. In particular, the ability of the modeling platform to represent the dynamics (temporal evolution) on a yearly or better monthly or even daily time scale is not considered in this study, while this would be essential to assess whether the modeling platform is able to meet its primary aim (i.e., the aforementioned motivations).
• section 4.1.
  • esp. l. 347-348: Do the authors use the longest available period for each gauge or the longest overlapping period (i.e., max 9 years between 2002 and 2010)? In any case, this relies on the hypothesis that an observed average over a few years (sometimes even only two years) as well as an observed average over two to several years covering another period (2002-2021) is representative for a steady state based on 1981-2010. I am not convinced that this hypothesis is true. I could agree that, the longer the observation period is, the closer the average gets to a steady state over the same period, even if this should still be verified. But in my opinion, there is no guarantee that the average over 2002-2021 is representative for the 1981-2010 steady state as this ignores potential shifts in the terrestrial water regime, e.g., due to climate change. One could think of the impact on streamflow of earlier snow melt, less snow accumulation in winter, more frequent extreme events, multi-year droughts, etc. Assuming that an average over a few years is representative for a steady state might be even more questionable, since these few years could be characterized by extreme events (droughts, floods). The resulting average would certainly not correspond to a steady-state, preventing a robust comparison. The authors already mention the potential impact of hydraulic engineering (e.g., dam operations) on the average streamflow (see l. 354), especially over short time periods.
• section 4.2.
  • the structure of this section might be improved to make it easier for the reader to follow. For example, after the first paragraph (l. 364-371) describing Figure 6, one would expect a first analysis of these results. But this analysis only starts on l. 432.
  • l. 364 and Figure 6: Are the steady states of the two global datasets over the same period as for CONCN (i.e., 1981-2010)? If not, is the hypothesis valid that these steady states, which might have been reached under different climatic conditions, are comparable? For example, if one region experiences less (or more) precipitation and/or higher evapotranspiration due to climate change, the resulting steady state will very likely be different.
  • l. 397-416: In the same way as my comments above for the evaluation of streamflow, I do not see any reason why one could assume that the observed average over 2018 could be considered as representative or close to a steady state generated with data from 1981-2010. Especially for water table depth with a potentially huge impact of inherited conditions from previous years (memory effect), not only 2018, but also the previous years would need to be close to the 1981-2010 average hydrologic regime to – maybe – approach a steady-state-like state. I understand that this is the reason why the authors try to strengthen their evaluation with the analysis of the residuals in the context of the long-term trend from GRACE, thereby trying to make a link between the steady state based on 1981-2010 and 2018, but, in my opinion, the uncertainty of this whole evaluation remains high. The only way to provide a robust and representative evaluation is to do the comparison over a common period.
  • l. 401-404: This might be easier to understand if the authors could briefly explain why this analysis integrating GRACE data is needed.
  • l. 411: if the global models are calibrated, do they not implicitly account for human interaction, via the observational data used for calibration? This would then be contradictory with the statement in l. 406.
  • l. 421-422: In l. 377, I understand that these grid cells are excluded from this analysis (precisely for the reason explained here). Please clarify this.

Specific comments

• l. 54-55: While I agree that it is pressing to develop a modeling platform accounting for it, this statement suggests that CONCN accounts for water quality control, which is not the case.
• l. 99: About the “unique dramatic topographic relief”. On one side, each part of the world has a unique relief, thus I could agree with this formulation. On the other side, many other regions (e.g., the US, South America, Africa, Europe, New Zealand, Japan, etc.) have transitions from mountains to coastal plains, thus facing similar challenges for hydrologic modeling.
• Figure 1: What is the meaning of the white coloring within the model domain in Fig. 1f? Here, I would interpret it as “no data”, is that correct? If it is zero, it should be colored according to the color bar (i.e., dark blue). If it is “no data”, how do the authors deal with it as source-sink term for ParFlow? This should be clarified in the text.
• l. 167-169: The procedure is not clear to me. Did the authors generate D8 connectivity slopes in addition to the aforementioned D4 slopes? If yes, why was it needed? What do they mean by “vector networks”?
• l. 180: How are the sinks handled? Is the inflowing ponding water removed before/after each time step?
• l. 200: Why do the authors derive the soil texture from this global dataset instead of using directly the soil hydraulic properties?
• l. 208: What is meant by “flow barriers”? Is the permeability in these grid cells further reduced by a factor or set to a very small value?
• l. 250: “several locations” It might be useful to add some more details here: How many locations? Are they distributed over the country and/or hydroclimatic regions to ensure a representative analysis and selection of datasets?
• l. 259: Why did the authors use the period 1981-2010 and not, e.g., 1991-2020? This could have made the evaluation easier, as the authors state further below that more observations are available for the last years (esp. since the 2000s).
• l. 275: It might be useful to add a source or some additional explanation on how Manning’s coefficient is “set to vary by land cover type”.
• l. 284: What is the role of the seepage face boundary condition?
• l. 284-285 and 287: On which time scale does the total storage change have to be less than 1% (resp. 3%)? Is it e.g., between two consecutive time steps or on an inter-annual basis?
• l. 292: I understand that it is important to reach an equilibrium for groundwater and for river discharge, but does a quasi-steady state for discharge in arid and semi-arid regions really make sense? Is the resulting discharge not too far away from reality? I would guess that in reality, the discharge is highly variable in these regions, with very low flow, or even no flow at all, most of the time alternating with high discharge after precipitation events or snow melt.
• Figure 3: How can the authors explain that they have streamflow values everywhere and not just in the streambeds in the south-east and north-east of the model domain? Or is it just an impression due to the visualization of a dense hydrographic network?
• Figure 3: It would be useful to add in the caption which period is shown. Or is it the end of the spinup (i.e., resulting quasi-steady state)?
• Figure 5: It might be useful to add in the caption that the gauges are grouped per basin as shown on Fig. 1b.
• Figure 7 and in the text: Do the “residuals” correspond to the difference between CONCN and the observed values at the wells?
• l. 453: All regions in the world experience increasing extreme weather events such as droughts and floods. What may make China “one of the most significant ecohydrologic hotspots in the world” could be the intense water use in the highly populated areas of the country. However, this is not accounted for in the model platform presented in this paper.

Technical corrections

• l. 29: Meaning of RSR?
• l. 56: Correct “ with a 10 km resolution”.
• l. 67: Meaning of USGS?
• l. 112: Correct “key components of the ParFlow model”?
• Figure 1: The north-eastern edge of the domain is hidden behind the color bars.
• Figure 1: In the caption, what do “f.g.”, “sil.”, and “c.g.” stand for?
• l. 154-156: For clarity, it might be good to specify that this concerns each grid cell individually, e.g., something like “D4 connectivity means that, within each grid cell, streamflow is allowed…”.
• l. 182: Correct “with those in IHU”?
• l. 199: In l. 125, the thickness of the second layer (from the top) is 0.3 m. Here, it is indicated to be 0.4 m.
• l. 260: Correct “Tarim River Basin”?
• l. 267: It is important to expand the acronym of CLM to avoid any confusion, as nowadays CLM usually means “Community Land Model”, while the CLM integrated in ParFlow is the “Common Land Model”.
• l. 340-343: There is a mismatch in the number of gauges: 95 (total) – 6 (no location) – 1 (close to another) – 1 (outside of domain) = 87, not 88.
• Figure 7: Correct “The background shows the average decrease of groundwater storage”.
• l. 397: Correct “by the three models”?
• Figure 8: Indicate in the caption that you compare the steady state over 1981-2010 with observations from 2018.
• l. 433: Correct “and the two global models”?
• l. 435: Correct “across the three models”?
• l. 445: Correct “below – these require” or maybe “below. These require”?
• l. 515: Correct “have been cited”?
• l. 524: Correct “reported in this paper”?
• l. 525: Correct “which is a consortium”?
• l. 526: Correct “and the Office”?

Yang et al. present results for a coupled surface-groundwater model for continental China. The work is appropriate for HESS, the research is very interesting and of high quality, and the manuscript is well written. Nonetheless, I have several recommendations regarding the evaluation of the model(s). I believe that a more process-oriented evaluation would be more meaningful for both the authors and the readers. My main other concern right now is a lack of discussion of the results. Both aspects are straightforward to rectify though.

Larger Comments:

[1] The use of a scaled statistical error metric: The authors state that “Note that all performance evaluations in this paper are based on the RSR value which is the ratio of the root mean squared error to the standard deviation of observations. An RSR value of 1.0 suggests good performance while 0.5 suggests excellent performance (O'neill et al., 2021).” These qualitative statements go back to the paper by Moriasi et al. (2007, doi.org/10.13031/2013.23153) who suggested some subjective qualification for normalized statistical metrics. The use of this subjective language persists even though it has been shown multiple times that the ease with which such values can be achieved varies with system properties (e.g. DOI: 10.1002/hyp.6825; doi.org/10.5194/hess-23-4323-2019). Therefore, these statements of good or poor performance with fixed thresholds are very unhelpful because – depending on the system modelled – it will be easy or hard to achieve these values. Personally (the authors do not have to share this view), I find it much more helpful to assess which system properties allow for high or low model performances (e.g. DOI 10.1088/1748-9326/abfac4 Figure 3 or DOI 10.1088/1748-9326/ad52b0). Such analyses are particularly valuable when done across multiple models, which often show that many models work well under specific conditions (often high wetness levels).

[2] Possibility for understanding process controls: The focus on statistical metrics and maps for the comparison of the model with observations or other models provides limited insights into how and (potentially) why the models differ. A simple but effective way to provide more insight is to plot the water table depth (WTD, or other output variables) against (potentially) controlling variables as functional relationships. For example, when plotting WTD against topographic slope for two of the models used by the authors – GLOBGM and Fan, the recent study by Reinecke et al. (doi.org/10.1088/1748-9326/ad8587) showed that GLOBGM is strongly correlated with slope, while the Fan model and global observations do so much less. Also, the Fan model shows distinct WTD differences between water and energy limited regions, while GLOBGM hardly does so. Similarly to my point 1, what controls the variability of model outputs and the output differences? These plots would include data, which the authors should have readily available – hence there is not much additional effort needed to try this.

[3] Model omissions: Over 0.5 million km2 of Southern China has Karst geology (doi.org/10.1007/s10980-019-00912-w), which shows significantly different recharge patterns than many other geologies (doi.org/10.1073/pnas.1614941114). How is this reflected in the model set-up? Do these regions show distinctly different patterns than other areas regarding recharge or other variables?

[4] Comparison with global models: Global models are rather crude approximations of local hydrology – shown regularly. Comparison to these models is a good starting point, but also limited in what one can learn. Do any national scale modelling efforts exist for China that would also provide a comparison for the model introduced here? Clearly the model presented here has tremendous potential – given its coupled nature – but how would it have to be further improved? It would be interesting to discuss more what additional aspects local or regional models might consider relevant.

[5] Lack of discussion: As is often the danger when Results and Discussion sections are not separated, there is a lack of actual discussion. The discussion section should place the results in context of existing literature. This has not yet been done. Other evaluations of the models used exist. Other modelling studies have assessed different strategies for China or globally Etc. The authors need to place their results into such context, preferably by separating Results and Discussion into distinct sections.

Minor Comments:

[6] Line 85ff.: The authors state that “Significant progresses or consensus have been achieved in community discussions regarding model parameterization, evaluation, calibration, and intercomparison”. Given that at least the cited Gleeson et al. stresses the current lack of adequate evaluation strategies for global models, I would personally not frame it quite this positively. I do think that there is still significant advancement needed to derive at adequate strategies, and I also think that consensus is not yet there.

[7] Figure 6. The lower plots show positive and negative deviations from 0. The maps would be much clearer if the authors were to use a diverging color scheme as they do in Figure 7.  Though I can also see that the authors prefer to keep the colors similar to the actual values.

Thorsten Wagener