Natural neighbour co‐ordinates (Sibson co‐ordinates) is a well‐known interpolation scheme for multivariate data fitting and smoothing. The numerical implementation of natural neighbour co‐ordinates in a Galerkin method is known as the natural element method (NEM). In the natural element method, natural neighbour co‐ordinates are used to construct the trial and test functions. Recent studies on NEM have shown that natural neighbour co‐ordinates, which are based on the Voronoi tessellation of a set of nodes, are an appealing choice to construct meshless interpolants for the solution of partial differential equations. In Belikov et al. (Computational Mathematics and Mathematical Physics 1997; 37(1) : 9–15), a new interpolation scheme (non‐Sibsonian interpolation) based on natural neighbours was proposed. In the present paper, the non‐Sibsonian interpolation scheme is reviewed and its performance in a Galerkin method for the solution of elliptic partial differential equations that arise in linear elasticity is studied. A methodology to couple finite elements to NEM is also described. Two significant advantages of the non‐Sibson interpolant over the Sibson interpolant are revealed and numerically verified: the computational efficiency of the non‐Sibson algorithm in 2‐dimensions, which is expected to carry over to 3‐dimensions, and the ability to exactly impose essential boundary conditions on the boundaries of convex and non‐convex domains. Copyright © 2001 John Wiley & Sons, Ltd.