Review Lung et al. (2023) - Open Boundary Conditions for Atmospheric Large Eddy

Simulations and the Implementation in DALES4.4

The manuscript describes the numerical implementation of inflow and outflow boundary conditions into the LES model DALES. By using an idealized big-brother like LES with periodic boundary conditions for a convective boundary layer, the performance of the newly implemented boundary conditions was systematically evaluated.

The topic itself fits well into the journal and the papers contains some interesting thoughts about the implementation of the in- and outflow boundary conditions as well as its evaluation.

However, I have major concerns about the simulation setup, the conceptualization how the method is evaluated, as well as the description of the method.

Hence, I cannot recommend the manuscript for publication until major revisions have been made. My major concerns and minor comments are outlined in the following.

Major comments

A - The author motivate their study by nesting LES domains into large-scale model domains. It is well known that other LES models which use Dirichlet boundary conditions for time-dependent mesoscale flow inputs sometimes suffer from wave-like structures near the boudaries, so better formulated boundary conditions to overcome this would be highly appreciated. However, as far as I understand, the boundary conditions described herein are only supposed to be used for idealized situations where the inflow and outflow boundary are fixed over the LES simulation period. For example, in a mesoscale-nested simulation, it is likely that the wind speed and direction continuously change in time, meaning that an inflow boundary can become an outflow boundary and so on. While this is still considered in the equations, though not supported by any analysis, the situation where a lateral boundary can become both, inflow and outflow boundary at the same time, is not considered in the equations. For example, this situation can occur if you want to model mesoscale phenomena like sea breezes, local wind systems, convective situations with weak winds, or situations like frontal passages. This is because the radiation boundary condition requires slab averages of the outward-pointing component. If there is a significant inflow at this boundary, the <u_n> can become negative. In case this happens, the flow becomes quickly unstable in conjunction with radiation boundary conditions, meaning that the proposed method is only applicable for idealized scenarios. Thus, the use of a slab average actually prohibits that a boundary can be both, inflow and outflow boundary at the same time. I recommend to rephrase the general motivation in this context, in order to avoid the impression that the proposed formulation of the boundary conditions solves the issue in general.

B - The description of the boundary conditions lacks important information and is partly misleading. For example, the boundary conditions are formulated as tendencies instead of boundary values. However, the boundary value itself is required for the spatial descretization of the advection term, so I recommend to reformulate the equation towards boundary values. Further, the term slab average is not fully defined. It seems to have a different meaning at the outflow boundary compared to the inflow boundary. Moreover, the formula for the time-scale computation seems to be wrong because the second term in Eq. 13 does not become dimensionless.

C - The setup description of the test case lacks important information. Which surface boundary conditions did the authors use (momentum, heat, SGS-TKE, ...), which numerical schemes were applied (pressure solver, advection and time discretization, ...). Moreover, it is not clear to me how the north and south boundaries were treated (period BC vs. inflow/outflow BC?).

D - I like the idea of a big-brother simulation to investigate the impact of the open boundary conditions in a systematic manner. However, the performance of the open BC is not sufficiently supported by the test case and the analysis. The authors only used a single setup for a convective boundary layer with a fixed inflow and outflow boundary. However, convection may easily masked systematic effects because instantaneous fluctuations may superimpose weaker systematic biases. For this purpose I think the evaluation of the model need to be extended towards purely neutral flows. Moreover, I think the test scenario should be also extended to a case with changing inflow conditions with respect to the wind speed to i) evaluate the performance of the mass-conservation scheme and ii) to demonstrate that proposed time-dependent relaxation time-scale algorithm works properly. Also a test case with changing wind direction is required to demonstrate that the boundary conditions can also deal with such situations.

Minor comments

l8: The first part of the sentence sounds strange and should be rephrased.

l12: I wouldn't say LES exists to study small scale weather phenomena but would formulate this in a more general way, e.g. to study turbulent motions.

l25: What do the authors mean by the term "fields"?

l43-45: It would be useful for the reader if the authors would be more specific, i.e. which model uses which kind of BC. The way the sentence if phrased is too general in my opinion. Also, concerning a description of inflow/outflow BC, the Maronga et al. (2015, https://doi.org/10.5194/gmd-8-2515-2015) paper is more suited reference.

l48-l49: In addition to the Mazzaro paper it would be nice to add the original literature (Mirocha et al., 2014, https://doi.org/10.1175/MWR-D-13-00064.1, plus the follow-up literature - see also references in Mazzaro et al., 2017) of the cell perturbation method too. Also, to my knowledge, Heinze et al. (2017) used no prescribed boundary conditions as stated in the follow sentence but periodic boundary conditions in combination with a large-scale forcing term inferred from mesoscale model output.

Intro: The manuscript would profit if the authors add some more text to introduce the term "open BC" and distinguish it from period boundary conditions with respect to its advantages and disadvantages. For example, also with periodic boundary conditions you can study larger-scale phenomena, even over heterogeneous land surfaces in particular cases.

l53: What do the authors mean with the term "numerical boundary layer"?

l54-55: The reference to Heinze et al. (2017) at this point is misleading and not correct. As mentioned before, the study used period BC and the relaxation therein does not refer of a relaxation in space but in time, formulated as a nudging term.

Moreover, a formulation like "often accompanied" is inappropriate here. The authors should be more specific in terms which model uses which strategy to mitigate boundary effects.

2.1.2 Inflow: What does it exactly mean that the Dirichlet condition is implemented as a tendency term? Suppose there is a mesoscale model input which changes over time and the LES model is in between 2 mesoscale model timesteps. How exactly are the BCs for the velocity vector and other quantities computed? I guess at the end DALES requires some kind of boundary values for each prognostic quantity for the spatial discretization rather than a tendency term?

Moreover, as the authors mentioned that a tendency term work well with the pressure solver, at what stage are the boundary values imposed, before or after invoking the pressure solver?

l170: I disagree with this interpretation. The boundary values enter the equations via the resolved- and subgrid-scale advection terms and not via the pressure term.

Headings of 2.1 and 2.2: The logical structure is misleading or the heading is poorly phrased. When 2.1 is about boundary-normal velocity components, I would expect that 2.2 is about boundary-parallel components and not about cell centered variables.

Equations - general: punctuation is missing

l190: What does the term "homogeneous Neumann condition" exactly mean? I see it is defined later in Eq. 11, but should be mentioned already when first used.

l200-201: To my knowledge this is exactly what is done in PALM (see Hellsten et al., 2021; Kadasch et al., 2021) and in WRF (Moeng et al., 2007; Mirocha et al., 2014), which does not seem to cause significant problems in both models. At least the authors should mention this. Furthermore, this raises the need to improve the argumentation why special Robin boundary conditions are required in conjunction to what happens in DALES when large gradients occur at the boundaries.

l227: Isn't e usually being defined as the SGS-TKE? If yes, the units do not match (term in brackets needs to be dimensionless). If not, how is a subgrid-velocity being defined? SGS-models usually give estimations for the SGS-TKE but not for the velocities. There are formulations for SGS-velocites (see e.g. Weil et al., 2004; Weil, J.C.; Sullivan, P.P.; Moeng, C.H. The Use of Large-Eddy Simulations in Lagrangian Particle Dispersion Models. J. Atmos. Sci. 2004, 61, 2877–2887), but I have the impression that the authors mean something different.(?)

l244-245: Can the authors please specify if this is their personal experience, or if it is experience deduced from previous studies? To my knowledge, the current state of literature does not support to make such a statement - there exists no extensive quantitative comparison between different methods so far.

Also, I strongly doubt that temperature fluctuations give perse a better solution than just adding perturbations onto the velocity components because the physical mechanisms of turbulence development differ and might not fit to the physical setup. For instance, in purely-shear driven flows this can lead to long persisting streak-like structures.

l263-264: For demonstrating the benefit of a newly developed method it is inappropriate to say that other test cases are not shown because they yield similar good results. Either you have conducted these tests and show some results of them, or you don't. In my opinion, purely neutral tests give different insights in the performance of a method as just a convective case. Same with cloudy boundary layers, where it is not straightforward how cloud prognostic quantities provided by mesoscale scales are treated in the LES at the boundaries.

l265: Can the authors please be more specific? A w* = 1.5m/s can be achieved in different ways, e.g. by altering the surface flux or the boundary-layer depth. What was the prescribed heat flux in the simulations and how was the initial profile of potential temperature being defined?

l267: Do the authors have arguments why they used such an anisotropic grid?

l268: Was the dt really fixed to 5s? In a CBL the vertical component can become about 10 m/s. In conjunction with a dz = 20, time steps of 2s would be required to maintain numerical stability of the advection equation.

l270 and following: If I understand right, you did perform a forcing where the open BC LES is driven by a period LES. In this regard, it is not clear to me how the coupling was realized. Did you take spatially resolved data, or did you only took horizontal mean profiles? Did you prescribed boundary values at all lateral boundary, i.e. the east, west, north, south and top boundary, or only that the west boundary? I might be wrong, but according to Fig. 3 it looks like you used periodic BC along y. So my question: Does the north/south boundary act as inflow/outflow boundary at the same time? Does the left inflow boundary could be also an outflow boundary (in a CBL with 3m/s mean wind this can happen)? Same with the right "outflow" boundary.

I strongly recommend the authors revise the setup description and add more details to allow for a better understanding what was done.

Furthermore, I am interested how the authors realized the coupling technically (some note in the text might be nice). Was is realized by an offline approach where the data is stored in a separate file or via an MPI coupling strategy between the big-brother and the open-BC simulation?

Fig. 2: It would be easier to understand if you show absolute values rather than differences. Further, did you compute the profiles from the entire xy-domain or did you exclude some areas near the boundaries? In my opinion it does not make much sense to include areas where the flow is potentially affected by the boundaries because this can bias the result, even if the flow features in the interior of the model domain perfectly match.

l338-339: To thoroughly evaluate this, xz cross-sections are required. It could well be the case the authors just randomly picked a height which is only weakly affected, while other heights show significant up- or downdrafts near the boundaries.

l339-340 and Fig. 4: Resolved or subgrid TKE? In the first case, how did you calculate the TKE (formula, time-averaging of the total fluxes, etc.)? In the next sentence you mention that the TKE is averaged over half an hour, which partly answers my question, but I have the impression that the calculation of TKE is not completely correct in this case. According to what you wrote, you computed instantaneous values of TKE from sum_i < (u_i'(t))^2 > and average these over time. This only works when u' refers to a phase average where homogeneous conditions along y apply. However, if the north/south boundaries are also in/outflow boundaries, this is strictly speaking not the case. Alternatively, u_i'^2 can be computed via a time average.

caption Fig. 4: How can a black line indicate a "fixed" ratio? I guess you mean something like ratio between horizontal and vertical advection?

l352-354: Which data was exactly used for the wavelet analysis? Did you use a spatial or a temporal data series for the wavelet analysis. In the latter case, at which distance from the inflow boundary? Did you use timeseries at at single point of time dependent yz cross section data. Which mother wavelet was employed? More specific information is required.

355: I do not understand why the analysis window is outside the domain. Actually the hatched area is defined by the cone-of-influence in the wavelet literature, describing the area in the scalogramm which is not affected by boundary effects. The sentence should be rephrased accordingly.

l363 and following: I agree, but this is not surprising as you simply forced an LES with output from another LES under idealized conditions (no changing wind direction, not much change in mass flux, etc.). The authors should put their statements into the context what their test case really shows.

l392-397: This is an interesting point because it systematically investigates the overshooting of turbulence also seen in previous studies (Munosz-Esparza and Kosovic, 2018 - https://doi.org/10.1175/MWR-D-18-0077.1 ; Kadasch et al., 2021). I would encourage the authors to also discuss their findings in the context of previous studies.