**In a recent article in JAMES (Pressel et al., 2017), we explore how numerical error and subgrid-scale modeling in LES interact to determine the quality of LES of stratocumulus clouds and show that a technique called implicit large eddy simulation provides particularly high fidelity LES. Here we offer a bit of background and a discussion of that work. If you are not familiar with stratocumulus clouds you can see a high resolution LES of stratocumulus here.**

Stratocumulus clouds play an important role in modulating Earth’s climate. They are highly reflective to sunlight and preferentially occur over the oceans, where they significantly reduce the amount of solar radiation that would otherwise be absorbed by the underlying dark surface. So effective are stratocumulus in altering Earth’s energy budget that small changes in their aerial coverage could have impacts similar in magnitude to the impacts of all anthropogenic emissions of greenhouse gases (e.g. Slingo (1990) ). This makes understanding the dynamics of stratocumulus clouds and how they change with climate of pressing importance.

Yet making progress is challenging, in part because their dynamics involve three-dimensional, turbulent fluid motions on the scale of a few meters, the nonlinear interaction of turbulence with the formation of cloud droplets, and the interaction of these cloud droplets with radiation (e.g. Wood (2012)). For example, the formation of cloud droplets depends on the availability of sufficient water vapor for condensation to occur; without condensation, there is no cloud. Having enough water vapor depends upon sufficiently vigorous turbulent transport of water vapor from the ocean surface to the cloud layer, with the primary motive force of this turbulence coming from density currents driven by strong radiative cooling. The strength of radiative cooling in turn depends on the characteristics of the clouds formed by the condensed vapor. If any of these interrelated processes breaks down, the stratocumulus clouds disappear.

These strong interactions between small-scale and seemingly disparate physical processes make understanding the intricacies of stratocumulus dynamics from observations a nearly insurmountable challenge. Thus, research has increasingly focused on a technique called Large Eddy Simulation (LES) to numerically simulate stratocumulus clouds, based on the well understood equations that govern small-scale turbulent motions, moist thermodynamics, and radiative transfer. In our group, as discussed in a previous blog post, we have undertaken a significant effort to develop a new Open Source infrastructure, called PyCLES, for LES of a broad range of clouds types, including stratocumulus.

However, even using state-of-the-art techniques, LES of stratocumulus presents a challenge, in part because the nonlinearities present in real stratocumulus are also present in the LES. Solving the physical equations that govern stratocumulus inherently involves approximation. For example, the governing equations must be approximated, or discretized, so that they can be solved on a finite number of grid points called a mesh, which must be small enough to be stored in memory on a supercomputer. This approximation carries with it two points of difficulty:

- Some of the turbulent motions are so small that they cannot be captured on a mesh that can be stored in the memory of even the largest supercomputers, requiring that the effects of the unresolved motions must be approximated.
- The process of writing the equations in such a way that they can be solved on a discrete mesh carries with it numerical error.

These approximation errors interact with the nonlinearities inherent in stratocumulus dynamics and can severely degrade even the bulk statistics predicted by LES.

The difficulty that LES has in simulating stratocumulus is apparent in the Stevens et al. (2005) LES intercomparison, which sought to explore the robustness of LES in simulating stratocumulus by comparing the results of a large number of LES to high quality observations. The results were striking: across the intercomparison the cloud fraction was found to vary by a factor of almost five and the amount of liquid water in the clouds by nearly a factor of twelve, with the simulations having the highest cloud fraction and liquid water agreeing most closely with the observations.

In our recent paper, Pressel et al. (2017), we sought to uncover the role played by the way in which the governing equations are solved on the mesh and the models used to represent unresolved turbulent processes in determining the quality of LES. To do this, we ran a large set of simulations of the Stevens et al. (2005) case in which we varied the numerical methods used to approximate the governing equations, as well as the models used to represent unresolved turbulent processes. PyCLES, the code developed in our group, includes a diverse set of numerical methods which makes it well equipped for this type of study.

What we find is that nearly the entire intercomparison range seen in Stevens et al. (2005) could be achieved in PyCLES by varying only the numerical methods and the models for unresolved turbulence. By classifying the simulations based upon the properties of the numerical methods used and the approach to representing unresolved turbulence we could begin to understand how these interact to affect the fidelity of the LES.

This classification is based upon how the numerical errors in the methods used to solve transport equations for momentum and scalar quantities (like entropy and specific humidity) affect small-scale gradients of these fields. In particular, PyCLES includes two types of transport schemes: central schemes, which tend to lead to overprediction of small-scale gradients, and weighted essentially non-oscillatory (WENO) schemes, which tend to dissipate small-scale gradients. The overprediction of small-scale gradients which occurs in central schemes tends to generate oscillatory numerical solutions. WENO schemes avoid this oscillatory behavior but introduce numerical dissipation. As the most widely used models for unresolved turbulence dissipate small scale gradients, the dissipation inherent in WENO schemes allows PyCLES to be run with no explicit model for unresolved processes, an approach that is commonly referred to as implicit large eddy simulation (ILES).

LES-predicted cloud fractions are shown below, along with results from the Stevens et al. (2005) intercomparison. The observed cloud fraction for this case was 1.0.

**Figure 1:Time series of cloud fraction for each of the experiments. The absolute and interquartile ranges of the Stevens et al. (2005). LES intercomparison are shown in light and dark shading, respectively. The solid back line shows the Stevens et al. (2005) intercomparison mean. (Figure and caption are from Pressel et al. (2017)) .**

The numbers in the figure legends refer to the accuracy of the numerical schemes. What might at first be quite surprising is that the type of numerical scheme is generally much more important to maintaining physically realistic values of cloud fraction than the accuracy of the scheme. In the figure, central schemes are indicated by even numbers (2, 4, 6) while WENO schemes correspond to odd numbers (5, 7, 9). The first number in each grouping gives the scheme used for solving the momentum equations and the second number identifies the scheme used for solving the scalar equations. The inclusion of models for unresolved turbulent motions is indicated by “M” and/or “S”, for momentum and scalar, respectively, while “N” means that no model for unresolved turbulence is used. Note that the “N” option, which corresponds to the implicit LES approach, is numerically stable only in conjunction with WENO schemes, whereas central schemes must be paired with a turbulence model.

A few points stand out:

- Mixing central schemes for momentum and WENO schemes for scalars (plots in the left-hand column) gives generally poor performance, although results are substantially improved when the WENO schemes are employed in the implicit LES approach.
- Pairing schemes (plots in the right-hand column) gives cloud fractions at or near one, except for the 22MS case, which appears to be degraded by the low numerical accuracy of the second-order scheme.
- The best results (Paired-NSGS plot) are obtained using a fully implicit approach to LES in which both momentum and scalar equations are solved using WENO schemes and no model for unresolved turbulence is employed.

Thus careful consideration of how LES are configured is necessary to produce high fidelity simulations of stratocumulus. Surprisingly, we find that using WENO numerics in an ILES (Paired-NSGS) configuration provides optimal results, given the configurations tested. This means that the dissipative numerical error inherent in the WENO scheme serves as an effective model for unresolved turbulent processes.

In Pressel et al. (2017) we conclude that the effectiveness of WENO schemes comes about through their ability to reduce mixing between the cloudy boundary layer and the overlying, much drier, free troposphere. This reduced mixing with the free troposphere prevents the cloud layer from drying, keeping the amount of liquid water concentrations in the cloud layer high, cloud-top cooling intense, and vertical turbulent mixing between the surface and cloud layer vigorous. Thus the stratocumulus clouds are maintained. Moreover, we contend that the quality of the simulations using the Mixed-SGS approach are severely degraded by excessive dissipation that occurs when the dissipation of the scalar transport schemes is combined with the dissipation due to models for unresolved turbulence that depend on spuriously large velocity gradients produced by the central difference numerics used for momentum. This induces anomalous mixing between the cloudy boundary layer and the free troposphere, which thins the stratocumulus clouds. We also provide a simple analytical model that supports these arguments. For a deeper discussion of these results, we encourage you to read our paper.