AbstractThe scaled boundary finite element method is a novel semi‐analytical technique, whose versatility, accuracy and efficiency are not only equal to, but potentially better than the finite element method and the boundary element method for certain problems. This paper investigates the possibility of using higher‐order polynomial functions for the shape functions. Two techniques for generating the higher‐order shape functions are investigated. In the first, the spectral element approach is used with Lagrange interpolation functions. In the second, hierarchical polynomial shape functions are employed to add new degrees of freedom into the domain without changing the existing ones, as in the p‐version of the finite element method. To check the accuracy of the proposed procedures, a plane strain problem for which an exact solution is available is employed. A more complex example involving three scaled boundary subdomains is also addressed. The rates of convergence of these examples under p‐refinement are compared with the corresponding rates of convergence achieved when uniform h‐refinement is used, allowing direct comparison of the computational cost of the two approaches. The results show that it is advantageous to use higher‐order elements, and that higher rates of convergence can be obtained using p‐refinement instead of h‐refinement. Copyright © 2005 John Wiley & Sons, Ltd.