The accuracy, efficiency, and stability of three numerical models with application to open ocean problems

Abstract The inviscid barotropic vorticity equation is integrated under a variety of assumed initial and boundary conditions corresponding to linear and nonlinear box modes, forced nonlinear box modes, and linear and nonlinear Rossby waves (with and without mean advection). The former two classes of problems are defined within a closed domain. The latter is partially or totally open to a specified external environment and therefore represents prototype limited-are calculations for the ocean. To determine the extent to which the accuracy and efficiency of limited-area calculations depend on the numerical integration scheme, each test problem is solved independently using the finite-difference (FD), finite-element (FE), and pseudospectral (PS) techniques. The three numerical models differ primarily in the formal accuracy of their spatial approximations and their treatment of vorticity at outflow points along the boundary. The FD model employs a centered second-order differencing scheme and requires an extrapolatory (computational) boundary condition to fix the values of vorticity at outflow boundary points. The FE model, which represents ψ and ζ as a summation of piecewise linear elements, is of fourth order for the linearized one-dimensional advective equation. Further, a technique is developed by which the determination of the interior values of ζ is decoupled from that of the boundary values; hence, the vorticity boundary conditions can be implemented without iterative techniques. Lastly, the “infinite-order” PS model avoids the assumption of lateral periodicity by expanding ψ and ζ in a double series of Chebyshev polynomials. The resulting vorticity equation is solved in the spectral domain using a modified alternating direction implicit method. All three models are of second order in time and have conservative formulations of the nonlinear terms. Integrations of moderate length (5–10 periods of the known analytic solution) are performed to determine the accuracy, stability, and efficiency of each model as a function of problem class and the associated physical and computational nondimensional parameters. The most important of these parameters are ϵ, the Rossby number; ν, the number of spatial degrees of freedom (grid points, expansion functions, etc.) per half wavelength of the reference solution; and η, the number of time steps per period of the reference solution. The latter two parameters are nondimensional measures of the spatial and temporal resolution of the numerical approximation. These tests show that all three models are, in general, capable of delivering stable and efficient solutions to linear and weakly nonlinear problems in open domains (0 ⩽ ϵ ⩽ 0.4, 4 ⩽ ν ⩽ 10, 64 ⩽ η ⩽ 128). Despite their added complexity, however, the FE and PS models are on the average, 4 and 15 times more accurate, respectively, than the FD model even taking into account its increased efficiency. The results also suggest that given a judicious selection of a frictional (filtering) mechanism and/or computational boundary condition (to suppress the accumulation of grid-scale features), each of the models can be made similarly accurate for highly nonlinear calculations ( ϵ ⪢ 0.4).