Response to the reviewer comments Title : Comparison of dealiasing schemes in large-eddy simulation of neutrally-stratified atmospheric boundary-layer type flows

Reviewer general comment: The paper compares various approaches for the dealiasing of the non-linear terms in large eddy simulation of atmospheric flows. Given that spectral methods are widely used in studying such flows under idealized conditions since they offer higher speed and better accuracy, the general theme is of interest to GMD readers. The paper does a good job in presenting the fundamentals of the problem, the proposed solutions, and how they compare when implemented in an actual code. But major revisions are needed. In particular, the study is a valuable comparison of the methods that is not available (to the best of my knowledge) in the literature and the authors should not try to conclude that one is more optimal than the others. They can simply present their findings and let the users determine which method is suitable for their needs.

Based on the comments from both reviewers, we have added a new discussion section to improve the readability of the manuscript. We have also improved the results section by running additional simulations at the resolutions of 192, 192, 128 and 192 3 , as well as running the 256 3 -case long enough to be able to report converged statistics. As a result Figures 7, 8 and 9 have been updated, and a new figure ( Figure 10) has been included, which is discussed at the end of the results section.

Changes in manuscript:
The new discussion section is reported here: "In the development of this manuscript, focus has been directed to the study of the advantages and disadvantages of different dealiasing methods. In this regard, throughout the analysis we have tried to keep the structure of the LES configuration as simple and canonical as possible, to remove the effect of other add-on complexities. Additional complications might arise when considering additional physics; here we discuss the potential effect that these different dealiasing methods could have on them. One of such elements of added complexity is for example the use of more sophisticated subgrid scale models based on dynamic approaches to determine the values of the Smagorinsky constant (Germano et al., 1991;Bou-Zeid et al., 2005). In most of these advanced subgrid models, information from the small-scale turbulent eddies is used to determine the evolution of the subgrid constant. However, in both the FT and FS method, the small turbulent scales are severely affected and hence use of dynamic subgrid models could be severely hampered unless these are accordingly modified and adjusted, maybe via filtering at larger scales than the usual grid scale. Another element of added complexity consists in using more realistic atmospheric forcing, considering for example the effect of the Coriolis force with flow rotation as a function of height and velocity magnitude. In this case, we hypothesize that the FT method could lead to stronger influences on the resultant flow field as this dealiasing technique not-only affects the distribution of energy in the small turbulent scales, but also in the large scales (as apparent from Fig. 2), being these later ones potentially more affected by the Coriolis force. This represents a strong non-linear effect, that is hard to quantify and hence further testing, including realistic forcing with a geostrophic wind and Coriolis force would be required to better quantify these effects. Also often in LES studies of atmospheric flows one is interested in including an accurate representation of scalar transport (passive/active). In this case the differential equations don't include a pressure term and hence most of the computational cost is linked to the evaluation of the convective term. As a result, the benefit of using alternative, cheaper dealiasing techniques (FT or FS) will be even more advantageous, yet the total gain is not trivial to evaluate a priori, and the effect on the scalar fields should also be further evaluated.
In general, we believe that it is not fair to advocate for one or other dealiasing method based on the results of this analysis. Note that the goal of this work is to provide an objective analysis of the advantages and limitations that the different methods provide, letting the readers the ultimate responsibility to choose the option that will adjust better to their application. For example, while having exact dealiasing (3/2-rule) might be better in studies focusing on turbulence and dispersion, one might be well-off using a simpler and faster dealiasing scheme to run the traditionally expensive warm-up runs, or to evaluate surface drag in flow over urban and vegetation canopies, where most of the surface force is due to pressure differences (Patton et al.,2016)."

Specific responses
Major comments: 1. Reviewer comment: While the FS method seem to be giving an acceptable performance as the authors argue, I wonder whether the ABL LES community should be going in a direction of saving computing time rather than maximizing the accuracy of the computation. We push for higher resolution to gain better accuracy and, with increasing computing power, I wonder whether a 20% drop in simulation time is worth it. We use dynamic SGS models that increase the computing time by 20% all the time. The plots in Fig 7 do not indicate that the FS method is as good as the 3/2 method. So in general I think the authors should not focus on the conclusion that the FS method is a good surrogate. They should present the information and findings, which will help modelers decide on the trade offs they want (on my end this convinces me that using a 3/2 method is indeed worth it.).
Authors response: Thanks, we indeed agree with the reviewer's point. We believe that gaining a 20% in computational time can be of interest in certain occasions, for example during warm up periods. However, as the reviewer mentions the strength of this manuscript should reside on conveying the facts of using different dealiasing methods, and allowing the corresponding end users to decide what's best for them according to their application.

Changes in manuscript:
To clarify this point we have added a discussion section (reported above) and rewritten the conclusion section. The new conclusion section emphasize the trade offs of each method. The modifications to conclusion reported here: "The Fourier-based pseudo-spectral collocation method (Orszag, 1970;Orszag and Pao, 1975;Canuto et al., 2006) remains the preferred "work-horse" in simulations of wall-bounded flows over horizontally periodic regular domains, which is often used in conjunction with a centered finite-difference or Chebychev polynomial expansions in the vertical direction (Shah and Bou-Zeid,2014;Moeng and Sullivan, 2015). This approach is often used because of the high-order accuracy and the intrinsic efficiency of the fast-Fourier-transform algorithm (Cooley and Tukey, 1965;Frigo and Johnson, 2005). In this technique, aliasing that arises when evaluating the quadratic non-linear term in the NS equations can severely deteriorate the quality of the solution and hence need to be treated adequately. In this work a performance/cost analysis has been developed for three well-accepted dealiasing techniques (3/2-rule, FT and FS) to evaluate the corresponding advantages and limitations. The 3/2-rule requires a computationally expensive padding and truncation operation, while the FT and FS methods provide an approximate dealiasing by low-pass filtering the signal over the available wavenumbers, which comes at a reduced cost.
The presented results show compelling evidence of the benefits of these methods as well as some of their drawbacks. The advantage of using the FT or the FS approximate dealiasing methods is their reduced computational cost (∼15% for the 128 3 case, ∼25% for the 256 3 case), with an increased gain as the numerical resolution is increased. Regarding the flow statistics, results illustrate that both, the FT and the FS methods, yield less accurate results when compared to those obtained with the traditional 3/2-rule, as one could expect.
Specifically, results illustrate that both the FT and FS methods over-dissipate the turbulent motions in the near wall region, yielding an overall higher mass flux when compared to the reference one (3/2-rule). Regarding the variances, results illustrate modest errors in the surface-layer, with local departures in general below 10% of the reference value across the considered resolutions. The observed departures in terms of mass flux and velocity variances tend to reduce with increasing resolution. Analysis of the streamwise velocity spectra has also shown that the FT method redistributes unevenly the energy across the available wavenumbers, leading to an over-estimation of the energy of some scales by up to 100%. Contrary, the FS methods redistributes the energy evenly, yielding a modest +13% energy magnitude throughout the available wavenumbers. Compared to the 3/2-rule, these differences in flow statistics are the result of the sharp low-pass filter applied in the FT method and the smooth filter that characterizes the FS method." 2. Reviewer comment: How do the FT and FS method influence the potential use of dynamic models that require good accuracy on the smallest resolved scales? If as the spectra show they damp these scales, than that would preclude using dynamic models and would be a significant disadvantage of FT and FS. The authors have in their code some dynamic models, they could perform the dynamic computations while still using the Static Smagorinsky (compute a dynamic Cs but dont use it).
Authors response: Thank you, this is a great point. In this regard, we have now run some additional simulations using the dynamic Smagorinsky model. Results are illustrated in Figure  1. In this figure, it can be observed that the dynamic model fails to compute the Smagorinsky constant when using either the FS or the FT methods as traditionally implemented. In this case the FT method strongly suppresses turbulence, resulting in the laminarization of the flow. Alternatively, while the consequences of using the FS method are less dramatic, the flow also exhibits a large acceleration at the top of the domain. As mentioned by the reviewer, these results are not surprising given that the dynamic models are using a relation between the small scales to compute the Smagorinsky constant. We believe that these could though be slightly improved by using information from scales larger than the traditional filtering scale, if the reader was really interested. Yet this reamins outside the scope of this manuscript. Therefore, it is advisable not to use the FT and FS methods with dynamic SGS models. In this regard, we have added additional text in the new discussion section that relates to the use of dynamic models and the fact that they probably require some modification to run with the FT and FS methods.

Changes in manuscript:
A discussion related to this comment has been added to the discussion section (see new discussion and conclusion section in comment #1).
3. Reviewer comment: Fig 8d and the associate sentence "Interestingly, results of the vertical flux (or stress, resolved and SGS) of stream-wise momentum (figure 8(d)) illustrate a good agreement between the different scenarios." The authors should be careful in this interpretation. The constant pressure gradient forcing requires and forces the stress profile to be linear. Regardless of how turbulence ends up looking like the turbulent fluxes have to adjust to balance the mean ∂P/∂x. What this figure indicates is that the SGS fraction is not strongly affected by the choice of dealiasing method, which is a good thing.
Authors response: Thank you for bringing this to our attention. As mentioned above the results section has been adjusted where this is taken care of.
4. Reviewer comment: Figure 9, and more generally: I would have liked to see a direct comparison of the largest scales (by filtering all simulations at n∆ , where n correspond to start of the damping or cutoff in figure 1) to see if the differences are only on the smallest scales or not (although given the mean velocity profiles, I suspect they are not).

Authors response: Thank you very much, this is a very interesting and important point.
In order to answer this comment and clarify the effect of the FT and FS methods on the spectra, we have developed an additional analysis using the spectra presented in the paper. Although the effect of the FT and FS methods on the small scales can be clearly observed on the spectra, their effect on the large scales cannot be directly assessed from the figure, as pointed out by the reviewer. To compute a direct comparison scale by scale, the following ratio was used for the 128 3 , 192 3 , and 256 3 -simulations, where E u,k denotes the power spectral density of the u velocity component at wavenumber k, XX stands for the dealiasing method FT or FS. Hence, if ρ(k) < 0 energy is removed at that scale, and if ρ(k) > 0 energy is added at that scale. Figure 2 presents the ratio ρ(k) for both methods where it can be observed that the effect of FT methods is very large at all scales. The large scales (0 ≤ k/k max ≤ 0.2) are affected with a reduction of energy of ∼25%. The mid-range scales (0.2 ≤ k/k max ≤ 0.6), corresponding to the inertial subrange, exhibit an overestimation of their energy of about ∼50% on average. Therefore, this method redistributes the energy of the small scales into the inertial sub-range scales. On the contrary, in the FS method, the energy from the filtered small-scales is redistributed more or less uniformly throughout with an averaged overall variation of less the 13%.

Changes in manuscript:
We have added figure 2 and its interpretation in the manuscript. This now reads as: "Although the effect of the FT and FS methods on the small scale can be clearly observed in figure 9, their effect on the large scales also needs to be quantified. To compute a direct comparison scale by scale, the following ratio was used (equation 2) for the 128 3 , 192 3 , and 256 3 -simulations, where E u,k denotes the power spectral density of the u velocity component at wavenumber k, XX stands for the dealiasing method FT or FS. If ρ(k) < 0 the energy density at that given wavenumber (k) is less than the corresponding one for the run using the 3/2 rule, viceversa if ρ(k) > 0. Figure 2 presents the ratio ρ(k) for both methods.
When using the FT method, energy at the low wavenumbers is underpredicted, whereas energy at the large wavenumbers is overpredicted. Departures are in general larger with decreasing resolution, with an excess of up to 100% for the 128 3 -simulations in the wavenumber range close to the cutoff wavenumber. On the contrary, when using the FS method, the energy from the filtered (dealiased) small-scales is redistributed quasi-uniformly throughout the spectra with an averaged overall variation of less than 13%." Minor comments: 1. Reviewer comment: Title is long and too descriptive: how about replacing the wordy "atmospheric boundary-layer type flows" with "atmospheric flows". One in fact could foresee using such methods for cloud resolving LES outside the ABL. Same on last line of abstract: why restrict the applications only to ABL flows?
Authors response: Thank you for this comment. We have adapted these very interesting suggestions. The new title is included below, as well as the new abstract.

Changes in manuscript:
Title: "Comparison of dealiasing schemes in large-eddy simulation of neutrally-stratified atmospheric flows" Abstract: "Aliasing errors arise in the multiplication of partial sums, such as those encountered when numerically solving the Navier-Stokes equations, and can be detrimental to the accuracy of a numerical solution. In this work, a performance/cost analysis is proposed for widely-used dealiasing schemes in large-eddy simulation, focusing on a neutrallystratified, pressure-driven atmospheric boundary-layer flow. Specifically, the exact 3/2 rule, the Fourier truncation method, and a high order Fourier smoothing method are inter-compared.
Tests are performed within a newly developed mixed pseudo-spectral collocation -finite differences large-eddy simulation code, parallelized using a two-dimensional pencil decomposition. A series of simulations are performed at varying resolution and key flow statistics are inter-compared among the considered runs and dealiasing schemes. Both the Fourier Truncation and the Fourier Smoothing method correctly predict basic statistics. However, they both prove to provide less accurate flow statistics when compared to the traditional 3/2-rule. The accuracy of the methods is dependent of the resolution. The biggest advantage of both of these methods against the exact 3/2-rule is a notable reduction in computational cost with an overall reduction of 15% for a resolution of 128 3 , 17% for 192 3 and 21% for 256 3 ." miss-understanding.

Changes in manuscript:
The text now reads as: "Note that the molecular viscous term has been neglected within the flow. However, the effect of the molecular viscosity at the surface is modeled using the logarithmic law, where the surface drag is parameterized through the surface roughness." 15. Reviewer comment: Page 6 line 25: what does "module" mean? Do they mean modulus?
Authors response: Indeed. We have clarified this in the manuscript.
16. Reviewer comment: Should equation 12 include f i to be consistent ?
Authors response: Thanks, this has been corrected.
17. Reviewer comment: Page 8, lines 10-15: Authors should clarify this is with the baseline 3/2 dealiasing I presume. Also how does the parallelization method impact these numbers?
Authors response: Thank you for pointing this out, we used the 3/2 as a baseline. In addition, an other comment on this regard was also made by an other reviewer, and hence we have included some clarification in the manuscript.
Here is our response to the question regarding the parallelization and the pressure solver.
We have noticed that the cost breakdown for the resolution of the convective term and the Poisson solver is also influenced by the pencil decomposition.
When treating the convective term with the pencil decomposition, the communication cost increases with respect to the traditional slice parallelization. In this case, a total of nine transpositions are needed to compute the convective term, significantly increasing the computational cost.
Opposite, the Poisson solver becomes faster when using the pencil decomposition in comparison to the slice parallelization. Note that in the pseudo-spectral method the horizontal directions (x and y) are treated in Fourier space and only the vertical direction (z) remains in physical space, therefore each mode in k x and k y become independent of each other. In this case the system of equations, originally of size n x × n y × n z becomes n x × n y systems of n z equations, making each vertical line in the domain independent. The pencil decomposition can take full advantage of this fact, making the resolution of the Poisson equation faster. Specifically, once the domain is transposed in the Z-pencil (square pipe aligned with the z-coordinate), the process of solving each of the n x × n y systems does not require any communication, making it very efficient, and limiting its cost to the transposition between the different pencils.

Changes in manuscript:
We realize that this is an important detail that should be also mentioned in the manuscript. For this reason we have included a couple of lines in section 3.3.
The text now reads as: "In addition, it is important to note that the low computational cost of the Poisson solver is related to the the use of the pencil decomposition, which takes full advantage the pseudo-spectral approach. Specifically, the Z-pencil combines with the horizontal treatment of the derivatives to make the implementation of the solver very efficient." 18. Reviewer comment: Page 10, lines 1-2: The z 0 they impose is 1cm, which corresponds more to a grass field than to a sparse forest of to a farmland. I suggest they check Brutsaerts books rather than to Stull for z 0 .
Authors response: We realized that there was a mistake in the actual value of z 0 . This one is actually of 0.1m, which corresponds to a sparse forest according to Stull and Brutsaert. We added the reference to Brutsaert's book in the manuscript.
19. Reviewer comment: Figure 5 and other are difficult to read. Why not use colors for the online version (Color is free with EGU, no?) Authors response: Because we couldn't find any information in this regard on the publisher web page we decided to go on black & white. If the editor confirms that color figures are free of charge, then we would be happy to change them 20. Reviewer comment: Page 11, lines 9-11: 30% drop in the convective term cost is good but I would not say it is significant. It would only be equivalent to about 20% drop in total computing time (given Fig 2), which would only be equivalent to a 5% reduction in the resolution. So I would remove significantly on line 9.

Authors response:
We have removed the "significantly" in the text. Additional detail in this regard has been provided earlier in comment 17.
Authors response: Thank you for pointing out the repetition.
23. Reviewer comment: Page12, line 9: correct the misspelling of "stream-wise" Authors response: Thank you for pointing out the misspelling.
24. Reviewer comment: Page 12 line 25, and page 14 line 10: "differentiated" is an unclear word. Please remove and clarify the two sentences.
Authors response: Thank you for pointing out these two sentences. The discussion of the results have been changed and these sentences have been removed.
25. Reviewer comment: Page 14 lines 15-16. There wont be any dispersive stresses in their simulations over homogeneous terrain so why mention them?
Authors response: In order to avoid any confusion in the discussion, we remove the mention to dispersive stresses in the text