Most Earth system models today use one of two first-order methods to calculate and apply fluxes of conserved quantities (e.g. water, energy, and momentum) at the Earth’s surface. The explicit flux coupling method, calculates atmosphere-surface fluxes only once per coupling time step, based solely on the current model state. This method is easily implemented and allows model components a great deal of flexibility in how they apply those fluxes, but it can result in numerical instability when the coupling time step size is too large. Meanwhile, the implicit flux coupling method inserts part of the surface flux calculation between two phases of a solver for a band-diagonal system representing vertical diffusion. This method can be highly stable, but requires the atmosphere and surface calculations to use compatible time discretizations.

We evaluate alternative methods that allow flux coupling to be performed stably for longer time step sizes, but without requiring interleaving of the atmosphere PBL parameterization and surface energy balance calculations, or even for these parameterizations to be evaluated with the same time step size. These methods are applied to both a simplified model of the atmospheric boundary layer, and to the Energy Exascale Earth System Model version 3 (E3SMv3), focusing on momentum fluxes over land. We show that this allows a known numerical instability in E3SMv3 (caused by its use of explicit flux coupling) to be alleviated without major changes to timestepping strategy in either the atmosphere or land models.