the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
Numerically consistent budgets of potential temperature, momentum, and moisture in Cartesian coordinates: application to the WRF model
Matthias Göbel
Stefano Serafin
Mathias W. Rotach
Download
 Final revised paper (published on 26 Jan 2022)
 Preprint (discussion started on 20 Jul 2021)
Interactive discussion
Status: closed

RC1: 'Comment on gmd2021171', Anonymous Referee #1, 17 Aug 2021
This paper introduces a budget tool for the WRF model and shows a successful
closure of this budget with very small residuals based on a Cartesian
formulation of the advection terms. There are many details considered in
achieving this, especially how the timevarying heights of coordinate
surfaces in the masscoordinate model are represented. The paper is valuable
in describing this and WRFlux looks like it would be useful to the
community by saving rediscovering these details when others want to do
accurate budgets.The paper is very close to acceptable as is. I will note some minor issues
and some confusion regarding references to equation numbers that have to
be corrected for publication.Minor Points
1. Line 235. It would be good to mention the length of the simulation here as
this surface flux would only be reasonable for a period of a few hours.2. Figure 2 caption. Should this be Eq. 8 instead of 6?
3. Line 292295. Some orders of magnitude here seem ten times too large compared
to Figures 1 and 2. 10^(3) and 10^(4) look more accurate.4. Figure 4. Should be Eq. 13, 14 and 15.
5. Table 1. Also should be Eq. 13, 14 and 15.
6. Line 303. Is this the lowest level only?
7. Figure 5. Labels 15 and 16 should be 14 and 15. Caption appears correct.
8. Line 318. References to Eq. 15 and 14 should be 16 and 15.
9. Line 339. Maybe "total" is better than "final" here.
10. Figure 2. Comment: This large difference is interesting and I would like to
have a better conceptual idea of why. Is it because the Cartesian representation
is somehow less sensitive to coordinate motion? In Fig. 2b are we looking at the
expansion of the coordinate layers with heating?
AC1: 'Reply on RC1', Matthias Göbel, 20 Aug 2021
Thank you very much for your positive review and your comments!
Here are our replies:1. Line 235. It would be good to mention the length of the simulation here as this surface flux would only be reasonable for a period of a few hours.
We added the simulation time (4 h) in the following paragraph.
2. Figure 2 caption. Should this be Eq. 8 instead of 6?
Eq. 6 is correct. It refers to the budget equation in the terrainfollowing coordinate system as it appears in the WRF model equations.
3. Line 292295. Some orders of magnitude here seem ten times too large compared to Figures 1 and 2. 10^(3) and 10^(4) look more accurate.
The +/1e2 numbers refer to SGS and resolved contributions, respectively. These are not shown separately in any figure. Fig. 2 only shows the decomposition into total turbulent (resolved + SGS) and mean advection.
The total tendency close to the ground above the ridge can be seen in the solid green line in Fig. 3. It is about 0.04e3=4e5. This is why we wrote on the order of 1e5.
4. 5. and 7.:
Thanks for pointing us to the wrong equation references. After removing an equation, we forgot to update the hardcoded references in the figures and table. This is now fixed.
6. Line 303. Is this the lowest level only?
At the lowest mass level (~3m agl) the tendency in the green dashed line is positive, at the second level (~9m agl) it is strongly negative, and at the third level (~15m agl) it is strongly positive again. But also higher up the differences between the two green lines are rather large, at least up to the 15th vertical level (~100m agl). The AGL height on the yaxes of Fig. 3 and 4 had a small error in the computation. We have corrected these figures.
8. Line 318. References to Eq. 15 and 14 should be 16 and 15.
The references to Eq. 15 and 14 are correct here.
9. Line 339. Maybe "total" is better than "final" here.
ok
10. Figure 2. Comment: This large difference is interesting and I would like to have a better conceptual idea of why. Is it because the Cartesian representation is somehow less sensitive to coordinate motion? In Fig. 2b are we looking at the expansion of the coordinate layers with heating?
What we see in Fig. 2b is the LHS of Eq. 6. This is not the rho*theta tendency on constant model levels, but the tendency of rho*z_eta*theta (governing equation in WRF) normalized with mean(rho*z_eta).
If we would plot the rho*theta tendency on constant model levels (1st term on LHS of Eq. 8) it would look very similar to Fig. 2a. So yes, the difference between Fig 2a and 2b arises mainly from the expansion of the coordinate layers with heating (magnitude of z_eta increases with time). And yes, the representation in Fig. 2b is more sensitive to coordinate motion.
RC2: 'Reply on AC1', Anonymous Referee #1, 20 Aug 2021
Thanks for the responses and I am satisfied with them. I hope something can be added in the text about the difference in Fig. 2a and 2b to explain it as in this response, because this would be an advantage of using the Cartesian budget introduced to study physical tendencies that are separated from coordinate tendencies.

RC2: 'Reply on AC1', Anonymous Referee #1, 20 Aug 2021

AC1: 'Reply on RC1', Matthias Göbel, 20 Aug 2021

RC3: 'Comment on gmd2021171', Anonymous Referee #2, 03 Sep 2021

Review of "Numerically consistent budgets of energy, momentum and mass in Cartesian coordinates: Application to the WRF model" By Göbel et al. for GMD
This paper presents the design and implementation of an opensource software package called WRFlux for calculating the budgets of potential temperature, water vapor, momentum etc. on a Cartesian grid for WRF simulations. The paper is well written and is a valuable contribution to the modeling community that GMD serves. I recommend minor revision. Here are my comments:

As I understood it, the current implementation of WRFlux needs it to run online with WRF. I wonder if there is a way to provide an offline version without significant changes to the WRF code, which will be much easier to use for most people.

Section 2.1, Equation 11: I think it would be helpful to readers if a reference is made to the fact that all the WRF prognostic variables are socalled "coupled" (multiplied by the mass inside the grid cube per unit area) as explained in the technical notes?

Line 195: "The fluxes and all budget components except for advection are averaged in time during model integration". I don't know what "fluxes" refer to here. Surface fluxes?

Line 280: "The only difference between the total tendencies in the terrainfollowing and the Cartesian formulation is the second term on the lefthand side in Eq. 11, which accounts for the height of the vertical levels being timedependent". I wonder if the averaging contributes to the difference as well. For the first term in Eq. 11, the averaging is done after dividing by z_eta.

The equations numbers are probably wrong in the figure/table labels, shouldn't they be Eqs. 13, 14 and 15 rather than 14, 15 and 16?

Line 313: "To quantify the differences, we plot the sum of all forcing terms for each budget calculation method against the actual model tendency ..." By that sentence do you mean plotting the LHS and RHS of Eqs. 13, Eqs. 13 with 2ndorder advection, 14 and 15 respectively? If so, It would be clearer to say that explicitly.


AC2: 'Reply on RC3', Matthias Göbel, 17 Sep 2021
We are grateful to the reviewer for her/his insightful and constructive comments.
As we included a new equation and a new figure, the equation numbers and figure numbers have changed. Line, equation, and figure numbers in our replies refer to the revised manuscript.
 As I understood it, the current implementation of WRFlux needs it to run online with WRF. I wonder if there is a way to provide an offline version without significant changes to the WRF code, which will be much easier to use for most people.
Budget analysis on WRF simulations requires two steps: (1) Sending individual tendency terms to the WRF output stream; (2) Recursively averaging all relevant quantities over time to enable an estimation of the explicitly resolved turbulent fluxes. WRFflux implements the necessary code changes for both operations. Doing the same analysis entirely offline is not impossible in principle, but would not be efficient. As an alternative to step (1), one could develop code to compute budget terms from regular model output, mimicking exactly what WRF does; but this would imply rewriting large parts of WRF's dynamical core and parameterizations. As an alternative to step (2), one could output every single time step during the model integration and then do the averaging offline; but this would require immense storage space, besides being computationally inefficient in comparison to recursion. We believe our solution is an optimal compromise between the complexity of the task and the usability of the tool. The GitHub page (https://github.com/matzegoebel/WRFlux/) includes an extensive manual and we are open to help users getting the tool to run. 
Section 2.1, Equation 11: I think it would be helpful to readers if a reference is made to the fact that all the WRF prognostic variables are socalled "coupled" (multiplied by the mass inside the grid cube per unit area) as explained in the technical notes?
We added this in line 96.
 Line 195: "The fluxes and all budget components except for advection are averaged in time during model integration". I don't know what "fluxes" refer to here. Surface fluxes?
“Fluxes” refers to the subgridscale and resolved fluxes of the five prognostic model variables in the whole domain, not only at the surface. We changed this in line 212. The “except for advection” only refers to the budget components.
 Line 280: "The only difference between the total tendencies in the terrainfollowing and the Cartesian formulation is the second term on the lefthand side in Eq. 11, which accounts for the height of the vertical levels being timedependent". I wonder if the averaging contributes to the difference as well. For the first term in Eq. 11, the averaging is done after dividing by z_eta.
The coupled tendency is first calculated with the timeaveraged flux in Eq. 10 and then divided by <ρ z_{η}>. This was not stated clearly in the original manuscript, so we corrected the caption of Fig. 3.
 The equation numbers are probably wrong in the figure/table labels, shouldn't they be Eqs. 13, 14 and 15 rather than 14, 15 and 16?
Yes, indeed  thank you for spotting this. However, since we included an additional equation (Eq. 12) in the revised manuscript, the labels are correct now.  Line 313: "To quantify the differences, we plot the sum of all forcing terms for each budget calculation method against the actual model tendency ..." By that sentence do you mean plotting the LHS and RHS of Eqs. 13, Eqs. 13 with 2ndorder advection, 14 and 15 respectively? If so, It would be clearer to say that explicitly.
Yes. We modified lines 347348 and the caption of Fig. 6 accordingly.
 As I understood it, the current implementation of WRFlux needs it to run online with WRF. I wonder if there is a way to provide an offline version without significant changes to the WRF code, which will be much easier to use for most people.


RC4: 'Comment on gmd2021171', Anonymous Referee #3, 06 Sep 2021
This study addresses the important issue of budget analysis closure and provides a precise budget calculation tool that can benefit the numerical modeling community. When reading the abstract, one of my concerns is that the main concept and results (including the comparison with simplified approximations, i.e., using a lowerorder advection operator than the model) resemble those in Chen et al. (2020). However, after reading the entire manuscript, I think this article does contain substantial new content, and WRFlux has its unique usefulness and applicability. The practicability of Cartesian coordinate transformation is also clear in the example simulation with nonuniform orography and therefore strongly sloped η surfaces at low levels. In all, the advantage of using WRFlux is well exhibited, and the illustration materials are nicely presented. Still, I have some comments and suggestions:
Specific comments
 Title: Suggest changing "energy" and "mass" to “potential temperature” and “moisture”, respectively. Most results are about the θ tendency budget, which is related but not equivalent to “energy budget” in the strict sense. The current title may misleadingly imply an Earth’s energy budget, moist static energy budget, turbulence kinetic energy (TKE) or total energy budget, etc. Furthermore, I don’t think the result of the “mass budget” is ever shown.
 Abstract: I suggest reconstructing the abstract to emphasize the unique features and applicability of WRFlux to differentiate this work from Chen et al. (2020), such as the retrieval of the resolved turbulence component, transformation to the easytointerpret Cartesian coordinate, etc., all of which are particularly useful for largeeddy simulations or simulations over nonuniform topography.
 Introduction: It is important to mention that there have been some attempts to achieve a precise budget retrieval in the exact WRF model (Chen et al. 2020 is not the first), e.g., Lehner (2012, https://collections.lib.utah.edu/ark:/87278/s61n8fxw), Moisseeva and Steyn (2014, https://doi.org/10.5194/acp14134712014) and Potter et al. (2018, https://doi.org/10.1029/2018JD029427), etc.
 L3538: I’m not sure how the advective form may “hinder interpretation when internal processes dominate the budget of the control volume, which can be a single grid cell or a larger volume”? Can you please elaborate?
 Chen et al. (2020) noted that a closed budget for w is challenging, partially because w equation is implicit, coupled with the geopotential tendency equation, and partially due to the variable’s inherent rapid variation on small scales. How does WRFlux overcome this issue for the w budget? It is also surprising to see that WRFlux performs the best for w compared to all the other variables and that the accuracy (in terms of NRMS and r99; Table 1) is much higher for the momentum than for the thermodynamic variables (θ and q_{v}) in your case. Do you have any idea why? Is this result casedependent or relevant to the physical design of WRFlux?
 It would help readers follow Section 2 more easily if the connections between subsections are made more explicit and some terminologies are explained more precisely.
 For example, I find it confusing that “Eq. 11 is equivalent to Eq. 8” (L124128), but “Eq. 8 cannot be closed numerically in contrast to Eq.11” (L128). This confusion may be related to how you define the term “close” here. Suppose you are referring to “budget closure”. In that case, my understanding is that it depends on the accuracy of the conducted budget calculation, not the format of the chosen equation to which the budget analysis is applied. Do you mean that because Eq. 11 is numerically consistent with the equation used in the model, the lefthandside terms of Eq. 11 represent more accurately the actual tendency simulated by the model? And therefore, a budget analysis using Eq. 11 can achieve a higher degree of budget closure? Please clarify.
 Also, I suggest reminding readers at the beginning of Section 2.4 that Eq. 11 (i.e., Eq. 13 after discretization) will be used for the precise budget calculation. It would also be helpful to specify which terms in Eq. 11 can be decomposed into the mean advective and turbulent components (all righthandside terms excluding S or all righthandside terms excluding only some partial S, i.e., the subgrid physics?)  L113L114: “Numerical consistent” > “Numerical consistency”. Also, placing this sentence (“Numerically… after discretization.”) here seems odd. You are neither showing the discretization nor discussing the budget closure in the following lines yet.
 L130: It is still unclear to me why recalculating the vertical velocity is desired. If Eqs. 6 and 9 are also used in WRF, doesn’t it make more sense to use WRF’s prognostic to be numerically consistent with the model? Unless this is relevant to the issue I mentioned in comment #4?
 L191194: Does this part follow exactly Chen et al. (2020)?
 L195: “The fluxes and all… are averaged in time during model integration”. Averaged over the entire simulation? Or over a userspecified time window? If so, what are your suggested value and the physical reason behind it?
 L203204: ”Since WRF does not actually solve the continuity equation...” This may be misleading. WRF does solve the continuity equation, just not in the same format as expressed in your Eq. (12), i.e., in terms of the 3D variable ρ but column dry air mass per unit area (2D µ_{d}). Note that µ_{d} is indeed advanced with time in the WRF dynamical solver. I suggest changing the sentence to something like “Since WRF does not directly solve for ρ but the integrated column dry air mass…”
 L244: “The setup leads to … a very dry atmosphere, therefore moist processes are neglected.” Not sure if I overlook something but the initial setup of the moisture field is not mentioned?
 L246: “We calculate full θtendencies and decompose them into resolved turbulence, subgridscale turbulence and mean advective.” What about the rest of the retrieved budget terms, such as “physical parameterizations and numerical diffusion and damping”? To clarify, it may help to move the sentence “Since no microphysics scheme is activated and the simplified radiation scheme only affects the surface energy balance, the heat budget in the atmosphere only consists of resolved advection and subgridscale diffusion” (L261263) here and mention that for general applications, other gridresolved parameterized physics terms are possible and categorized as additional budget components.
 Figure 1: Panels ac show the total turbulence (trb), which is the sum of the resolved and subgridscale components. I’m interested in their individual contributions. E.g., What is the relative magnitude between these two components? Do they have similar spatial distributions with the same signs, or do they offset each other? Considering that one of the distinctive features of this tool from previous WRF budget retrieval works is the decomposition into mean and turbulent components, I suggest strengthening the relevant discussion.
 L283: “...in the alternative form of the equation (Eq. 8), the correction term for the time derivative is almost negligible” I’m lost. This is not shown or can be inferred by figure 2?
 L290295: Have you checked if the sum of all the budget components in this alternative analysis is in close balance with the sum of your Eq. 11(13) in Fig. 1?
Technical corrections:
 L113, L124, …: Replace “equation” with “Eq.” here. Please check the rest of the manuscript for consistency.
 L134: “energy” > “potential temperature”
 L154155: “Since… to derive.” This sentence is confusing. Suggest changing to “Although the momentum variables are staggered differently from the thermodynamic variables, their discretized equations can be derived analogously…”
 L226: “ridgetoridge” > ”bottomtobottom” or “valleytovalley” ?
 Figure 4 legend: I believe the legends “WRFlux (Eq. 14)” should be replaced by “WRFlux (Eq. 13)”, “Eq. 15” by “Eq. 14”, and “Eq. 16” by “Eq. 15”.
 Figure 5 subtitles for each panel: Same issue as above.
 Table 1: Same issue as above.

AC3: 'Reply on RC4', Matthias Göbel, 17 Sep 2021
We are grateful to the reviewer for her/his insightful and constructive comments.
As we included a new equation and a new figure, the equation numbers and figure numbers have changed. Line, equation, and figure numbers in our replies refer to the revised manuscript.
Specific comments

Title: Suggest changing "energy" and "mass" to “potential temperature” and “moisture”, respectively. Most results are about the θ tendency budget, which is related but not equivalent to “energy budget” in the strict sense. The current title may misleadingly imply an Earth’s energy budget, moist static energy budget, turbulence kinetic energy (TKE) or total energy budget, etc. Furthermore, I don’t think the result of the “mass budget” is ever shown.
OK, we changed the title accordingly.

Abstract: I suggest reconstructing the abstract to emphasize the unique features and applicability of WRFlux to differentiate this work from Chen et al. (2020), such as the retrieval of the resolved turbulence component, transformation to the easytointerpret Cartesian coordinate, etc., all of which are particularly useful for largeeddy simulations or simulations over nonuniform topography.
Good point. To emphasize these features, we added a sentence and moved the last sentence to a more prominent place (L59).

Introduction: It is important to mention that there have been some attempts to achieve a precise budget retrieval in the exact WRF model (Chen et al. 2020 is not the first), e.g., Lehner (2012, https://collections.lib.utah.edu/ark:/87278/s61n8fxw), Moisseeva and Steyn (2014, https://doi.org/10.5194/acp14134712014) and Potter et al. (2018, https://doi.org/10.1029/2018JD029427), etc.
OK, we added these references in the introduction in L2829.

L3538: I’m not sure how the advective form may “hinder interpretation when internal processes dominate the budget of the control volume, which can be a single grid cell or a larger volume”? Can you please elaborate? We added a sentence there: “For instance, if the fluxes on two opposing sides of a grid box are equal, the flux form yields zero net tendency, while this is not necessarily the case for the advective form.” (L3940)
Lee et al. 2004 looked at the budget of a larger control volume and noted that this is due to “internal processes” in the control volume. We don’t want to explain this further in the paper, so we dropped the unclear mentioning of the “internal processes”.

Chen et al. (2020) noted that a closed budget for w is challenging, partially because w equation is implicit, coupled with the geopotential tendency equation, and partially due to the variable’s inherent rapid variation on small scales. How does WRFlux overcome this issue for the w budget? It is also surprising to see that WRFlux performs the best for w compared to all the other variables and that the accuracy (in terms of NRMS and r99; Table 1) is much higher for the momentum than for the thermodynamic variables (θ and q_{v}) in your case. Do you have any idea why? Is this result casedependent or relevant to the physical design of WRFlux? Chen et al. found that closing the wbudget offline is difficult or impossible due to the mentioned issues, but the online budget calculation works fine for them. They also achieved a higher degree of budget closure for w than for horizontal momentum. Our budget calculations are also done online considering all budget components for w, including the important acousticstep contributions.
In contrast to Chen et al. (2020), our budget calculations are timeaveraged. To check the effect of the time averaging on the budget closure, we did two small test simulations with and without time averaging. Without time averaging, the NRMSE scores in the terrainfollowing coordinate system become considerably better for q_{v} and w and slightly worse for θ, u, and v.
We did not investigate why the budget closure in the shown simulations is higher for momentum than for θ and q_{v}. However, in some of our test simulations, it was the other way around, so it is definitely casedependent (background wind speed, moisture content, physical parameterizations used, etc.).

It would help readers follow Section 2 more easily if the connections between subsections are made more explicit and some terminologies are explained more precisely.  For example, I find it confusing that “Eq. 11 is equivalent to Eq. 8” (L124128), but “Eq. 8 cannot be closed numerically in contrast to Eq.11” (L128). This confusion may be related to how you define the term “close” here. Suppose you are referring to “budget closure”. In that case, my understanding is that it depends on the accuracy of the conducted budget calculation, not the format of the chosen equation to which the budget analysis is applied. Do you mean that because Eq. 11 is numerically consistent with the equation used in the model, the lefthandside terms of Eq. 11 represent more accurately the actual tendency simulated by the model? And therefore, a budget analysis using Eq. 11 can achieve a higher degree of budget closure? Please clarify.
Eq. 11 is similar to WRF’s governing equations because the coordinate metric z_{η} appears within the derivatives. In the manuscript, we showed how the budget can be closed with Eq. 11. When using the product rule to transform Eq. 11 to Eq. 8, two terms cancel out analytically, but not numerically. This is because one of the canceling terms comes from a horizontal derivative and the other from a vertical derivative. Therefore, the righthand side of Eq. 16 is not identical to the righthand side of Eq. 14. Thus, we argue that Eq. 16 cannot be closed numerically in WRF. The lefthand sides of Eq. 14 and 16, instead, are almost identical in our simulations.
To clarify this point, we added Eq. 12 and a bit of explanation in L167170.
 Also, I suggest reminding readers at the beginning of Section 2.4 that Eq. 11 (i.e., Eq. 13 after discretization) will be used for the precise budget calculation. It would also be helpful to specify which terms in Eq. 11 can be decomposed into the mean advective and turbulent components (all righthandside terms excluding S or all righthandside terms excluding only some partial S, i.e., the subgrid physics?)
OK. We added a few introductory sentences in Section 2.4 with the suggested content (L184189).

L113L114: “Numerical consistent” > “Numerical consistency”. Also, placing this sentence (“Numerically… after discretization.”) here seems odd. You are neither showing the discretization nor discussing the budget closure in the following lines yet.
That’s right. Nevertheless, in our opinion, this sentence is important for the logical train of thought. We added a reference to the discretization section in L116.

L130: It is still unclear to me why recalculating the vertical velocity is desired. If Eqs. 6 and 9 are also used in WRF, doesn’t it make more sense to use WRF’s prognostic to be numerically consistent with the model? Unless this is relevant to the issue I mentioned in comment #4?
Recalculating w is necessary to obtain a closed budget in the Cartesian coordinate system. Originally, there was a significant difference in the calculation of z_{η}ω in WRF‘s implementation of Eq. 6 and 9. We reduced this difference strongly by changing the discretization of Eq. 9 in WRF (L216221, https://github.com/wrfmodel/WRF/pull/1338). Nevertheless, there is still a small difference remaining: In Eq. 6, z_{η} is included in the coupled prognostic variable µ_{d}ω, whereas in Eq. 9 it is calculated directly from the geopotential. To overcome this difference and obtain a closed budget, we recalculate w in a consistent way and use it in the vertical advection term. Whether this is appropriate is definitely debatable. But the effect on w itself is hardly noticeable. Only when adding up the large and counteracting budget components, the difference becomes noticeable. With the prognostic w, the budget closure (NRMSE) is about one order of magnitude worse (L356359). We suppose the reviewer wanted to ask whether this is relevant to the issue mentioned in comment #5, not #4. These issues are not connected. The recalculation of w mentioned is only relevant for the transformation to the Cartesian coordinate system. It is not used as the prognostic variable in the conservation equation for w. We clarified this in L217219.

L191194: Does this part follow exactly Chen et al. (2020)?
The retrieval of subgridscale tendencies is different since we only extract the fluxes and then do the tendency calculation offline (as for the resolved fluxes). The retrieval of other forcing terms (physics, pressure gradient, damping, diffusion, Coriolis) was implemented very similarly to Chen et al. but independently from them.

L195: “The fluxes and all… are averaged in time during model integration”. Averaged over the entire simulation? Or over a userspecified time window? If so, what are your suggested value and the physical reason behind it?
Over a userspecified time window (we added this in L213). In this study, we use 30 min (L272) as is often used for Reynolds averaging (compromise between large sample for statistics and stationarity). We added a sentence with a reference on this in L272274.

L203204: ”Since WRF does not actually solve the continuity equation...” This may be misleading. WRF does solve the continuity equation, just not in the same format as expressed in your Eq. (12), i.e., in terms of the 3D variable ρ but column dry air mass per unit area (2D µ_{d}). Note that µ_{d} is indeed advanced with time in the WRF dynamical solver. I suggest changing the sentence to something like “Since WRF does not directly solve for ρ but the integrated column dry air mass…”
Actually, WRF's mass continuity equation in terms of µ_{d} (Eq. 2.12 in the technical notes version 4, http://www2.mmm.ucar.edu/wrf/users/docs/technote/v4_technote.pdf) does correspond to the rightmost term in our Eq. 13, except for a factor g. In fact, µ_{d} and ρ_{d} are related by a diagnostic equation: µ_{d} = g ρ_{d} z_{η}. Ensuring that the last term of our Eq. 13 is exactly zero is difficult because the mass continuity equation in WRF is indeed solved but not integrated explicitly in the acoustic time steps, but rather µ_{d} is diagnosed from its definition (Eq. 3.20 in the technical notes) after vertically integrating the mass divergence (Eq. 3.22). We acknowledge that our original formulation was misleading, and we replaced it with the information given in the preceding sentence (L223225).

L244: “The setup leads to … a very dry atmosphere, therefore moist processes are neglected.” Not sure if I overlook something but the initial setup of the moisture field is not mentioned?
Yes. We added the specification of the moisture profile (constant RH=40%) in L263.

L246: “We calculate full θtendencies and decompose them into resolved turbulence, subgridscale turbulence and mean advective.” What about the rest of the retrieved budget terms, such as “physical parameterizations and numerical diffusion and damping”? To clarify, it may help to move the sentence “Since no microphysics scheme is activated and the simplified radiation scheme only affects the surface energy balance, the heat budget in the atmosphere only consists of resolved advection and subgridscale diffusion” (L261263) here and mention that for general applications, other gridresolved parameterized physics terms are possible and categorized as additional budget components.
We moved the sentence as suggested and included the suggested additional sentence in L268270.

Figure 1: Panels ac show the total turbulence (trb), which is the sum of the resolved and subgridscale components. I’m interested in their individual contributions. E.g., What is the relative magnitude between these two components? Do they have similar spatial distributions with the same signs, or do they offset each other? Considering that one of the distinctive features of this tool from previous WRF budget retrieval works is the decomposition into mean and turbulent components, I suggest strengthening the relevant discussion.
The subgridscale component is mainly relevant very close to the surface and has a similar magnitude but opposite sign as the resolved turbulence component. We added a profile plot showing this decomposition (new figure 2) and some explanation in L294300.

L283: “...in the alternative form of the equation (Eq. 8), the correction term for the time derivative is almost negligible” I’m lost. This is not shown or can be inferred by figure 2?
Yes, this is not shown (we added “not shown” in L316). The plot for Eq. 8 looks the same as Figure 3a. But the correction terms on the LHS of Eq. 8 and 11 are conceptually different (L134136) since the prognostic variable is different in both equations (ρθ vs ρz_{η}θ).

L290295: Have you checked if the sum of all the budget components in this alternative analysis is in close balance with the sum of your Eq. 11(13) in Fig. 1? Fig. 1 does not show the total tendency, only the turbulent and mean advective components, not the sum of both. The “total” in the plot refers to the sum of horizontal and vertical components. We added this clarification in the caption.
The total tendency is shown in Fig. 3a. Close to the surface on the ridge the tendency is about 0.03 x 10^{3} so on the order of 10^{5} K s^{1} as written in L329 and visible in Fig. 3a.
Technical corrections:

L113, L124, …: Replace “equation” with “Eq.” here. Please check the rest of the manuscript for consistency.
Thanks, we fixed this at several locations.

L134: “energy” > “potential temperature”
OK

L154155: “Since… to derive.” This sentence is confusing. Suggest changing to “Although the momentum variables are staggered differently from the thermodynamic variables, their discretized equations can be derived analogously…”
OK

L226: “ridgetoridge” > ”bottomtobottom” or “valleytovalley” ?
OK, we changed it to “valleytovalley”.

Figure 4 legend: I believe the legends “WRFlux (Eq. 14)” should be replaced by “WRFlux (Eq. 13)”, “Eq. 15” by “Eq. 14”, and “Eq. 16” by “Eq. 15”.
Yes. However, since we included an additional equation (Eq. 12) in the revised manuscript, the labels are correct now.

Figure 5 subtitles for each panel: Same issue as above.
Yes

Table 1: Same issue as above.
Yes
