How we compute
Exascale and performance
A method is only useful if it runs. We work on making high-order finite element codes fast on the machines that actually exist: many-core CPUs, GPUs, and the exascale systems being built now, where the hardware changes faster than scientific software can usually follow.
The largest supercomputers now reach a billion billion operations per second, but almost none of that capability is easy to use. Performance comes from vectorisation, from many-core parallelism and from GPUs, and each of those needs different programming models, hardware-specific optimisations, and data to be laid out differently. A code written for one generation of machine is routinely slower than it should be on the next.
High-order methods are unusually well suited to modern computing platforms. This is because the work per element is large and highly structured, and the operations at the heart of the method can be recast as dense linear algebra, which is exactly what modern hardware is built to do quickly. Realising that in practice is the research problem: matrix-free formulations that never build the global matrix at all, kernels vectorised across elements rather than within them, and implementations that reach a good fraction of peak on hardware.
The portability question sits alongside it. We have compared the main programming models for this (e.g. OpenMP, OpenACC and Kokkos) on the same problem rather than in isolation, because the interesting question is not whether a model works but what it costs in performance and in the readability of the code that scientists have to maintain for a decade. That concern also drives our work on automatic code generation, where the kernels are produced from a description of the method rather than written by hand for each architecture.
A third thread changes what is parallelised. Spatial parallelism eventually saturates: past a certain number of processors there are too few elements each for the communication to be worth it. Parallel-in-time methods find concurrency along the time axis instead, which is counter-intuitive for an initial value problem and is the point.
Much of this work has been performed under the ExCALIBUR programme, which is the UK’s effort to get its scientific software ready for exascale, but continues with many of our different funding bodies.
Accepted version, CC BY, fig. 5
Published version, CC BY, fig. 6
Funded by
-
RAEng Research Chair
A five-year chair developing the next generation of solvers for nuclear fusion modelling, following on from the NEPTUNE projects.
Read more -
NEPTUNE (2022)
Developing high-order solvers for nuclear fusion modelling under the strategic priorities fund ExCALIBUR programme.
-
Exa-UQ
Uncertainty quantification for exascale computing.
-
APinTA
Parallel-in-time methods for exascale computing.
-
Addressing mesh generation and adaptation as a high-priority exascale use case under phase 1 of ExCALIBUR.
Read more
Who works on this
All of the team
Prof. David Moxey
Group leader and project leader of the Nektar++ framework.
Previously
Dr. Jacques Xing
Parallel-in-time methods for high-order methods.
Now: Research Software Engineer, Department of Aeronautics, Imperial College London
Recent papers
All publications-
Enhancing the Nektar++ spectral/hp element framework for parallel-in-time simulations
Comput. Phys. Commun., vol. 312, pp. 109584
Nektar++ is an open-source spectral/hp element framework designed to support the development of solvers for partial differential equations (PDEs). The software supports various discretization techniques, notably continuous, discontinuous, and hybridizable discontinuous Galerkin projections, in combination with both modal and nodal expansions. The framework uses object-oriented programming paradigms and modern C++ to encapsulate the complexity of the discretization techniques and contains several pre-implemented solvers, while new solvers can be easily implemented in customized proxy-apps. Nektar++ has been designed for high-scalability on high-performance computing (HPC) clusters using distributed memory architecture using a message passing interface (MPI). Recently, time-parallel integration techniques are being recognized as a potential solution to further increase concurrency and computational speed-up beyond the limits of strong scaling obtained from a pure spatial domain decomposition. Amongst the various time-parallel approaches proposed in the literature, the Parareal algorithm is a non-intrusive and iterative approach, exploiting a fine and a coarse solvers to achieve time-parallelism, and can be applied to both linear and non-linear problems. The efficient implementation of the Parareal algorithm in the Nektar++ open-source framework is described in this work, including a discussion of specific techniques used to adapt the code to a time-parallel framework. We demonstrate the application of these methods to multiple linear and non-linear problems provided by the existing Nektar++ solvers.@article{xing-2025, title = {Enhancing the Nektar++ spectral/hp element framework for parallel-in-time simulations}, author = {Xing, J. and Cantwell, C. D. and Moxey, D.}, journal = {Comput. Phys. Commun.}, year = {2025}, volume = {312}, doi = {10.1016/j.cpc.2025.109584}, pages = {109584}, url = {https://www.sciencedirect.com/science/article/pii/S0010465525000876} } -
Comput. Phys. Commun., vol. 284, pp. 108624
We develop efficient kernels for elemental operators of matrix-free solvers of the Helmholtz equation, which are the core operations for incompressible Navier-Stokes solvers, for use on graphics-processing units (GPUs). Our primary concern in this work is the extension of matrix-free routines to efficiently evaluate this elliptic operator on regular and curvilinear triangular elements in a tensor-product manner. We investigate two types of efficient CUDA kernels for a range of polynomial orders and thus varying arithmetic intensities: the first maps each elemental operation to a CUDA-thread for a completely vectorised kernel, whilst the second maps each element to a CUDA-block for nested parallelism. Our results show that the first option is beneficial for elements with low polynomial order, whereas the second option is beneficial for elements of higher order. The crossover point between these two schemes for the hardware used in this study lies at around P=4-5, depending on element type. For both options, we highlight the importance of the layout of data structures, which necessitates the development of interleaved elemental data for vectorised kernels, and analyse the effect of selecting different memory spaces on the GPU. As the considered kernels are foremost memory-bandwidth bound, we develop kernels for curved elements that trade memory bandwidth against additional arithmetic operations, and demonstrate improved throughput in selected cases. We further compare our optimised CUDA kernels against optimised OpenACC kernels, to contrast the performance between a native and a portable programming model for GPUs.@article{eichstadt-2023, title = {Efficient vectorised kernels for unstructured high-order finite element fluid solvers on GPU architectures in two dimensions}, author = {Eichst\"adt, J. and Peir\'o, J. and Moxey, D.}, journal = {Comput. Phys. Commun.}, year = {2023}, volume = {284}, pages = {108624}, url = {https://www.sciencedirect.com/science/article/pii/S0010465522003435}, doi = {10.1016/j.cpc.2022.108624} } -
29th IAEA Fusion Energy Conference, proceedings
This article provides a broad overview of Project NEPTUNE, which aims to create code for simulating edge plasma physics at scales approaching the exascale, and is run by the UKAEA as part of the UK’s ExCALIBUR framework. It is a representative survey of scope and associated outputs, rather than a detailed exposition; the structure is such as to describe the core components of NEPTUNE (efforts toward simulation of plasma fluid and kinetic effects, and their synthesis into exascale-ready combined codes) and then give an indication of the wide range of other work performed under the project.@inproceedings{threlfall-2023, title = {Software for Fusion Reactor Design: ExCALIBUR Project NEPTUNE: Towards Exascale Plasma Edge Simulations}, author = {Threlfall, EJ and Akers, RJ and Arter, W and Barnes, M and Barton, M and Cantwell, C and Challenor, P and Cook, JWS and Coveney, PV and Dodwell, T and others}, booktitle = {29th IAEA Fusion Energy Conference, proceedings}, year = {2023}, url = {https://eprints.whiterose.ac.uk/204359/1/IAEA_CN_316_2132.pdf}, organization = {York} } -
Fast barycentric-based evaluation over spectral/hp elements
J. Sci. Comp., vol. 90, pp. 78
As the use of spectral/hp element methods, and high-order finite element methods in general, continues to spread, community efforts to create efficient, optimized algorithms associated with fundamental high-order operations have grown. Core tasks such as solution expansion evaluation at quadrature points, stiffness and mass matrix generation, and matrix assembly have received tremendous attention. With the expansion of the types of problems to which high-order methods are applied, and correspondingly the growth in types of numerical tasks accomplished through high-order methods, the number and types of these core operations broaden. This work focuses on solution expansion evaluation at arbitrary points within an element. This operation is core to many postprocessing applications such as evaluation of streamlines and pathlines, as well as to field projection techniques such as mortaring. We expand barycentric interpolation techniques developed on an interval to 2D (triangles and quadrilaterals) and 3D (tetrahedra, prisms, pyramids, and hexahedra) spectral/hp element methods. We provide efficient algorithms for their implementations, and demonstrate their effectiveness using the spectral/hp element library Nektar++.@article{laughton-2022, title = {Fast barycentric-based evaluation over spectral/$hp$ elements}, author = {Laughton, E. and Zala, V. and Narayan, A. and Kirby, R. M. and Moxey, D.}, journal = {J. Sci. Comp.}, year = {2022}, volume = {90}, pages = {78}, url = {https://link.springer.com/content/pdf/10.1007/s10915-021-01750-2.pdf}, doi = {10.1007/s10915-021-01750-2} } -
8th European Congress on Computational Methods in Applied Sciences and Engineering
@inproceedings{liu-2022, author = {Liu, B. and Cantwell, C. D. and Moxey, D. and Green, M. and Sherwin, S. J.}, title = {Vectorised spectral/hp element matrix-free operator for anisotropic heat transport in tokamak edge plasma}, booktitle = {8th European Congress on Computational Methods in Applied Sciences and Engineering}, doi = {10.23967/eccomas.2022.291}, url = {https://www.scipedia.com/public/Liu_et_al_2022b}, year = {2022} }
Collaborators
- Spencer Sherwin Imperial College London
- Chris Cantwell Imperial College London
- Mike Kirby University of Utah
Interested in working on this?
We are always glad to hear from prospective PhD students and postdocs, and from groups who would like to work with us. Our openings page lists what is currently funded and the routes that are open year round.
Join us