Convergence testing is a common practice in the development of dynamical cores of atmospheric models but is not as often exercised for the parameterization of subgrid physics. An earlier study revealed that the stratiform cloud parameterizations in several predecessors of the Energy Exascale Earth System Model (E3SM) showed strong time step sensitivity and slower-than-expected convergence when the model's time step was systematically refined. In this work, a simplified atmosphere model is configured that consists of the spectral-element dynamical core of the E3SM atmosphere model coupled with a large-scale condensation parameterization based on commonly used assumptions. This simplified model also resembles E3SM and its predecessors in the numerical implementation of process coupling and shows poor time step convergence in short ensemble tests. We present a formal error analysis to reveal the expected time step convergence rate and the conditions for obtaining such convergence. Numerical experiments are conducted to investigate the root causes of convergence problems. We show that revisions in the process coupling and closure assumption help to improve convergence in short simulations using the simplified model; the same revisions applied to a full atmosphere model lead to significant changes in the simulated long-term climate. This work demonstrates that causes of convergence issues in atmospheric simulations can be understood by combining analyses from physical and mathematical perspectives. Addressing convergence issues can help to obtain a discrete model that is more consistent with the intended representation of the physical phenomena. Plain Language Summary Computer codes that simulate the time evolution of a physical system produce errors in the results due to the finite step sizes used to advance the calculations in time. These errors are expected to decrease as the time steps are shortened, at a rate determined by the characteristics of the equations and numerical methods. An earlier study revealed that the error reduction in several predecessors of the Energy Exascale Earth System Model (E3SM) was at a rate slower than expected. This study creates a simplified configuration of those models and investigates the causes of the unexpected behavior. We show that slow error reduction can be understood and improved by combining analyses from physical and mathematical perspectives. The required code modifications can lead to significant changes in the simulated long-term behavior of a full-fledged climate model. Furthermore, ensuring a proper rate of error reduction can help to obtain a computer code that is more consistent with the intended representation of the corresponding physical system.

Anomalous diffusion can be characterized by a mean-squared displacement x 2 (t) that is proportional to t a where a = 1. A class of one-dimensional moving boundary problems is investigated that involves one or more regions governed by anomalous diffusion, specifically subdiffusion (a < 1). A novel numerical method is developed to handle the moving interface as well as the singular history kernel of subdiffusion. Two moving boundary problems are solved: the first involves a subdiffusion region to the one side of an interface and a classical diffusion region to the other. The interface will display non-monotone behaviour. The subdiffusion region will always initially advance until a given time, after which it will always recede. The second problem involves subdiffusion regions to both sides of an interface. The interface here also reverses direction after a given time, with the more subdiffusive region initially advancing and then receding.

When solving a moving interface problem, the interface can be tracked using a variety of methods. The level set method captures the interface as an isocontour of a scalar level set function. The method has many advantages, including the ability to express many geometric quantities, such as the interface curvature, as derivatives of the level set. However, the numerically constructed level set function may not be smooth enough to compute the required derivatives. Furthermore, the method overall is known for having a high sensitivity to numerical dissipation. The former of these two shortfalls is addressed by augmenting the traditional level set equations with an explicitly tracked interface curvature. The curvature is then updated alongside the level set through an additional advection equation. The latter shortfall is addressed by combining a new velocity extension, that better maintains the signed-distance property of the level set, with a reconstruct-evolve-average approach to advancing the advection equations. The new approach is shown to have less mass loss (numerical dissipation) and better accuracy than comparable level set approaches. Three scenarios are investigated: an interface moving according to an external velocity field, an interface moving according to the interface curvature (mean-curvature flow), and the air-water interface of a water drop moving according to the curvature-dependent fluid velocity (surface-tension driven flow). point on the interface. Thus, the gradient of the SDLS is neither too shallow to see large movements in interface positions from small perturbations in the level set value nor too steep to influence truncation errors when taking derivatives of the level set. In general, there is not a closed-form expression for the SDLS of an arbitrary interface. Thus, it is numerically constructed using a reinitialization process. There are two main variations of reinitialization: PDE-based reinitialization, developed by Peng et al. [21], and the fast marching method (FMM) reinitialization, developed by Sethian [26]. PDE-based reinitialization obtains the SDLS by solving a timedependent PDE, whose steady state solution is the SDLS. In theory, one can get arbitrarily close to the SDLS by evolving the PDE further and further in time. In practice, however, the combination of a smooth signum function and numerical errors from each iteration can cause the zero contour to move during the reinitialization process, leading to errors in the interface location.FMM reinitialization, on the other hand, computes the SDLS by marching signeddistance values out from the interface. This typically involves first using high-order interpolation near the interface to seed the values, then solving the Eikonal equation in a directional fashion. While this method is generally less computationally costly and less prone to movement of the zero contour than that of PDE-based reinitialization, the resulting SDLS is usually less smooth. This can be an issue, for example, if the interface velocity depends on the curva...

Global atmospheric models seek to capture physical phenomena across a wide range of time and length scales. For this to be a feasible task, the physical processes with time or length scales below that of a computational time step or grid cell size are simplified as one or more parameterizations. Inadvertent oversimplification can violate constraints or destroy relationships in the original physical system and consequently lead to unexpected and physically invalid behavior. An example of such a problem has been investigated in the work of Wan et al. (2020, https://doi.org/10.1029/2019MS001982). This work addresses the issues at a more fundamental level by revisiting the parameterization derivation. A derivation of an unaveraged condensation rate in the unaveraged equations, sometimes referred to as subgrid equations, provides a clear description and more accurate quantification of the condensation/evaporation processes associated with cloud growth/decay, while avoiding simplifications used in earlier studies. A subgrid reconstruction (SGR) methodology is used to connect the unaveraged condensation rate with the grid cell averaged equations solved by the global model. Analyses of the SGR method and the numerical results provide insights into root causes of inconsistent discrete formulations and nonphysical behavior. It is also shown that the SGR methodology provides a flexible framework for addressing such inconsistencies. This work serves as a demonstration that when nonphysical behavior in a parameterization of subgrid variability is avoided through rigorous mathematical derivation, the resulting formulation can exhibit both better numerical convergence properties and significant impact on long-term climate. Plain Language Summary Simulation of global climate is a difficult task even for the most powerful supercomputers available today. Instead of directly simulating every physical process in the atmosphere, average effects are typically considered to make global climate simulations feasible. As an example, it is currently not possible to track the entire lifespan of every cloud on the planet over tens or hundreds of years. Instead, the average cloud formation and evolution is modeled using a "parameterization." A parameterization typically makes assumptions about the underlying physical behaviors represented by the average values. This work highlights issues that arise when those assumptions are overly simple or inconsistent. It shows how a mathematically rigorous approach to parameterization development can help avoid, diagnose, and correct those issues. Specifically, an existing large-scale condensation (nonconvective cloud formation) parameterization is examined and improved upon, with the result both capturing more realistic physical behavior and being better positioned to make use of upcoming advances in computational power.