There are many advantages to formulating an ice-sheet model in terms of a variational principle. In particular, this applies to the specification of boundary conditions, which might otherwise be problematic to implement. Here we focus primarily on the frictional basal sliding boundary condition in a non-Newtonian Stokes model. This type of boundary condition is particularly difficult because it is heterogeneous, requiring both a Dirichlet (no-penetration) condition normal to the bed, and a Neumann (frictional sliding) condition tangential to the bed. In general, Neumann conditions correspond to natural boundary conditions in a variational principle; that is, they arise naturally in the variational formulation and thus need not be explicitly specified. While the same is not necessarily true of Dirichlet conditions, it is possible to enforce a no-penetration condition using Lagrange multipliers within the variational principle so that the Dirichlet condition becomes a natural boundary condition. Thus, in the case of ice sheets, all relevant boundary conditions may be incorporated in the variational functional, making them particularly easy to discretize. For the Stokes model, the resulting basal boundary condition is valid for arbitrary topographic slopes. Here we apply the same methodology to the Blatter– Pattyn higher-order approximate model, which is ordinarily limited to small basal slopes by the smallaspect-ratio approximation. We introduce a modification that improves on the accuracy of the standard Blatter–Pattyn model for all values of the basal slope, as we demonstrate in the slow sliding regime for which analytical results are available. The remaining error is due to the effects of the small-aspect-ratio approximation in the Blatter–Pattyn model.