The spectral/hp element method combines the geometric flexibility of classical h-type finite element techniques with the desirable resolution properties of spectral methods. In this approach a polynomial expansion of order P is applied to every elemental domain of a coarse finite element type mesh. These techniques have been applied in many fundamental studies of fluid mechanics [31] and more recently have gained greater popularity in the modelling of wave-based phenomena such as computational electromagnetics [14] and shallow water problems [4] - particularly when applied within a Discontinuous Galerkin formulation.

There are at least two major challenges which arise in developing an efficient implementation of a spectral/hp element discretisation:

- implementing the mathematical structure of the technique in a digestible, generic and coherent manner, and
- designing and implementing the numerical methods and data structures in a matter so that both high- and low-order discretisations can be efficiently applied.

In order to design algorithms which are efficient for both low- and high-order spectral/hp
discretisations, it is important clearly define what we mean with low- and high-order.
The spectral/hp element method can be considered as bridging the gap between
the high-order end of the traditional finite element method and low-order end of
conventional spectral methods. However, the concept of high- and low-order discretisations
can mean very different things to these different communities. For example, the
seminal works by Zienkiewicz & Taylor [36] and Hughes list examples of elemental
expansions only up to third or possibly fourth-order, implying that these orders are
considered to be high-order for the traditional h-type finite element community.
In contrast the text books of the spectral/hp element community typically show
examples of problems ranging from a low-order bound of minimally fourth-order up to
anything ranging from 10^{th}-order to 15^{th}-order polynomial expansions. On the other
end of the spectrum, practitioners of global (Fourier-based) spectral methods [12]
would probably consider a 16^{th}-order global expansion to be relatively low-order
approximation.

One could wonder whether these different definitions of low- and high-order are just inherent
to the tradition and lore of each of the communities or whether there are more practical
reasons for this distinct interpretation. Proponents of lower-order methods might highlight
that some problems of practical interest are so geometrically complex that one cannot
computationally afford to use high-order techniques on the massive meshes required to capture
the geometry. Alternatively, proponents of high-order methods highlight that if the problem of
interest can be captured on a computational domain at reasonable cost then using high-order
approximations for sufficiently *smooth* solutions will provide a higher accuracy for a given
computational cost. If one however probes even further it also becomes evident that the
different communities choose to implement their algorithms in different manners. For
example the standard h-type finite element community will typically uses techniques
such as sparse matrix storage formats (where only the non-zero entries of a global
matrix are stored) to represent a global operator. In contrast the spectral/hp element
community acknowledges that for higher polynomial expansions more closely coupled
computational work takes place at the individual elemental level and this leads to the use of
elemental operators rather than global matrix operators. In addition the global
spectral method community often make use of the tensor-product approximations
where products of one-dimensional rules for integration and differentiation can be
applied.