Partial differential equations (PDEs) are used with huge success to model phenomena across all scientific and engineering disciplines. However, across an equally wide swath, there exist situations in which PDEs fail to adequately model observed phenomena, or are not the best available model for that purpose. On the other hand, in many situations, nonlocal models that account for interaction occurring at a distance have been shown to more faithfully and effectively model observed phenomena that involve possible singularities and other anomalies. In this article we consider a generic nonlocal model, beginning with a short review of its definition, the properties of its solution, its mathematical analysis and of specific concrete examples. We then provide extensive discussions about numerical methods, including finite element, finite difference and spectral methods, for determining approximate solutions of the nonlocal models considered. In that discussion, we pay particular attention to a special class of nonlocal models that are the most widely studied in the literature, namely those involving fractional derivatives. The article ends with brief considerations of several modelling and algorithmic extensions, which serve to show the wide applicability of nonlocal modelling.

There is more to the computation of pseudospectra than the obvious algorithm of computing singular value decompositions on a grid and sending the results to a contour plotter. Other methods may be hundreds of times faster. The state of the art is reviewed, with emphasis on methods for dense matrices, and a Matlab code is given.

Microstructure is a feature of crystals with multiple symmetry-related energy-minimizing states. Continuum models have been developed explaining microstructure as the mixture of these symmetry-related states on a fine scale to minimize energy. This article is a review of numerical methods and the numerical analysis for the computation of crystalline microstructure.

Multigrid methods are the fastest known methods for the solution of the large systems of equations arising from the discretization of partial differential equations. For self-adjoint and coercive linear elliptic boundary value problems (with Laplace's equation and the equations of linear elasticity as two typical examples), the convergence theory reached a mature, if not its final state. The present article reviews old and new developments for this type of equation and describes the recent advances.

The objective of molecular dynamics computations is to infer macroscopic properties of matter from atomistic models via averages with respect to probability measures dictated by the principles of statistical physics. Obtaining accurate results requires efficient sampling of atomistic configurations, which are typically generated using very long trajectories of stochastic differential equations in high dimensions, such as Langevin dynamics and its overdamped limit. Depending on the quantities of interest at the macroscopic level, one may also be interested in dynamical properties computed from averages over paths of these dynamics.

This review describes how techniques from the analysis of partial differential equations can be used to devise good algorithms and to quantify their efficiency and accuracy. In particular, a crucial role is played by the study of the long-time behaviour of the solution to the Fokker–Planck equation associated with the stochastic dynamics.

Finite Difference (FD) methods approximate derivatives of a function by local arguments (such as du(x) / dx ≈ (u(x + h) − u(x − h))/2h, where h is a small grid spacing) – these methods are typically designed to be exact for polynomials of low orders. This approach is very reasonable: since the derivative is a local property of a function, it makes little sense (and is costly) to invoke many function values far away from the point of interest.

Essentially non-oscillatory (ENO) and weighted ENO (WENO) schemes were designed for solving hyperbolic and convection–diffusion equations with possibly discontinuous solutions or solutions with sharp gradient regions. The main idea of ENO and WENO schemes is actually an approximation procedure, aimed at achieving arbitrarily high-order accuracy in smooth regions and resolving shocks or other discontinuities sharply and in an essentially non-oscillatory fashion. Both finite volume and finite difference schemes have been designed using the ENO or WENO procedure, and these schemes are very popular in applications, most noticeably in computational fluid dynamics but also in other areas of computational physics and engineering. Since the main idea of the ENO and WENO schemes is an approximation procedure not directly related to partial differential equations (PDEs), ENO and WENO schemes also have non-PDE applications. In this paper we will survey the basic ideas behind ENO and WENO schemes, discuss their properties, and present examples of their applications to different types of PDEs as well as to non-PDE problems.

We survey general techniques and open problems in numerical linear algebra on parallel architectures. We first discuss basic principles of paralled processing, describing the costs of basic operations on parallel machines, including general principles for constructing efficient algorithms. We illustrate these principles using current architectures and software systems, and by showing how one would implement matrix multiplication. Then, we present direct and iterative algorithms for solving linear systems of equations, linear least squares problems, the symmetric eigenvalue problem, the nonsymmetric eigenvalue problem, and the singular value decomposition. We consider dense, band and sparse matrices.

In this work we present a brief introduction to hierarchical bases, and the important part they play in contemporary finite element calculations. In particular, we examine their role in a posteriori error estimation, and in the formulation of iterative methods for solving the large sparse sets of linear equations arising from finite element discretization.

The Lanczos and conjugate gradient algorithms were introduced more than five decades ago as tools for numerical computation of dominant eigenvalues of symmetric matrices and for solving linear algebraic systems with symmetric positive definite matrices, respectively. Because of their fundamental relationship with the theory of orthogonal polynomials and Gauss quadrature of the Riemann-Stieltjes integral, the Lanczos and conjugate gradient algorithms represent very interesting general mathematical objects, with highly nonlinear properties which can be conveniently translated from algebraic language into the language of mathematical analysis, and vice versa. The algorithms are also very interesting numerically, since their numerical behaviour can be explained by an elegant mathematical theory, and the interplay between analysis and algebra is useful there too.

Motivated by this view, the present contribution wishes to pay a tribute to those who have made an understanding of the Lanczos and conjugate gradient algorithms possible through their pioneering work, and to review recent solutions of several open problems that have also contributed to knowledge of the subject.

We consider time-dependent (linear and nonlinear) Schrödinger equations in a semiclassical scaling. These equations form a canonical class of (nonlinear) dispersive models whose solutions exhibit high-frequency oscillations. The design of efficient numerical methods which produce an accurate approximation of the solutions, or at least of the associated physical observables, is a formidable mathematical challenge. In this article we shall review the basic analytical methods for dealing with such equations, including WKB asymptotics, Wigner measure techniques and Gaussian beams. Moreover, we shall give an overview of the current state of the art of numerical methods (most of which are based on the described analytical techniques) for the Schrödinger equation in the semiclassical regime.

V. L. Rvachev called R-functions ‘logically charged functions’ because they encode complete logical information within the standard setting of real analysis. He invented them in the 1960s as a means for unifying logic, geometry, and analysis within a common computational framework – in an effort to develop a new computationally effective language for modelling and solving boundary value problems. Over the last forty years, R-functions have been accepted as a valuable tool in computer graphics, geometric modelling, computational physics, and in many areas of engineering design, analysis, and optimization. Yet, many elements of the theory of R-functions continue to be rediscovered in different application areas and special situations. The purpose of this survey is to expose the key ideas and concepts behind the theory of R-functions, explain the utility of R-functions in a broad range of applications, and to discuss selected algorithmic issues arising in connection with their use.

Recent developments in numerical methods for solving large differentiable nonlinear optimization problems are reviewed. State-of-the-art algorithms for solving unconstrained, bound-constrained, linearly constrained and non-linearly constrained problems are discussed. As well as important conceptual advances and theoretical aspects, emphasis is also placed on more practical issues, such as software availability.

One of the most difficult problems in the numerical solution of ordinary differential equations (ODEs) and in differential-algebraic equations (DAEs) is the development of methods for dealing with highly oscillatory systems. These types of systems arise, for example, in vehicle simulation when modelling the suspension system or tyres, in models for contact and impact, in flexible body simulation from vibrations in the structural model, in molecular dynamics, in orbital mechanics, and in circuit simulation. Standard numerical methods can require a huge number of time-steps to track the oscillations, and even with small stepsizes they can alter the dynamics, unless the method is chosen very carefully.

General linear methods, as multistage multivalue methods, are the natural generalizations of linear multistep and Runge-Kutta methods. This survey contains a discussion of the traditional methods and a motivation for the general linear type of generalization. The new methods are introduced in terms of their formulation and the basic properties of consistency, stability and convergence. The order of general linear methods has to be looked at from a new point of view and it is shown how to use an algebraic structure (equivalent to B-series) to express conditions for a given order. Linear and non-linear stability for the new methods brings the theories for the classical methods into a comprehensive formulation and known results are outlined. Recently a number of subfamilies have been introduced and some of these are considered in detail. This applies in particular to methods with the property known as ‘inherent Runge-Kutta stability’. These seem to have prospects of yielding useful and efficient methods, and some progress towards their practical implementation is outlined. Finally, the relationship between stability functions and order of methods is discussed in a setting wide enough to include general linear methods as well as multiderivative methods, such as Obreshkov methods. The classical barriers due to Ehle, Daniel-Moore and Dahlquist (second barrier) all fit into a common pattern and these are explored in a general setting.

A computational framework is presented for integrating the electrical, mechanical and biochemical functions of the heart. Finite element techniques are used to solve the large-deformation soft tissue mechanics using orthotropic constitutive laws based in the measured fibre-sheet structure of myocardial (heart muscle) tissue. The reaction-diffusion equations governing electrical current flow in the heart are solved on a grid of deforming material points which access systems of ODEs representing the cellular processes underlying the cardiac action potential. Navier-Stokes equations are solved for coronary blood flow in a system of branching blood vessels embedded in the deforming myocardium and the delivery of oxygen and metabolites is coupled to the energy-dependent cellular processes. The framework presented here for modelling coupled physical conservation laws at the tissue and organ levels is also appropriate for other organ systems in the body and we briefly discuss applications to the lungs and the musculo-skeletal system. The computational framework is also designed to reach down to subcellular processes, including signal transduction cascades and metabolic pathways as well as ion channel electrophysiology, and we discuss the development of ontologies and markup language standards that will help link the tissue and organ level models to the vast array of gene and protein data that are now available in web-accessible databases.

In this review we discuss bifurcation theory in a Banach space setting using the singularity theory developed by Golubitsky and Schaeffer to classify bifurcation points. The numerical analysis of bifurcation problems is discussed and the convergence theory for several important bifurcations is described for both projection and finite difference methods. These results are used to provide a convergence theory for the mixed finite element method applied to the steady incompressible Navier–Stokes equations. Numerical methods for the calculation of several common bifurcations are described and the performance of these methods is illustrated by application to several problems in fluid mechanics. A detailed description of the Taylor–Couette problem is given, and extensive numerical and experimental results are provided for comparison and discussion.

Neural networks (NNs) are the method of choice for building learning algorithms. They are now being investigated for other numerical tasks such as solving high-dimensional partial differential equations. Their popularity stems from their empirical success on several challenging learning problems (computer chess/Go, autonomous navigation, face recognition). However, most scholars agree that a convincing theoretical explanation for this success is still lacking. Since these applications revolve around approximating an unknown function from data observations, part of the answer must involve the ability of NNs to produce accurate approximations.

This article surveys the known approximation properties of the outputs of NNs with the aim of uncovering the properties that are not present in the more traditional methods of approximation used in numerical analysis, such as approximations using polynomials, wavelets, rational functions and splines. Comparisons are made with traditional approximation methods from the viewpoint of rate distortion, i.e. error versus the number of parameters used to create the approximant. Another major component in the analysis of numerical approximation is the computational time needed to construct the approximation, and this in turn is intimately connected with the stability of the approximation algorithm. So the stability of numerical approximation using NNs is a large part of the analysis put forward.

The survey, for the most part, is concerned with NNs using the popular ReLU activation function. In this case the outputs of the NNs are piecewise linear functions on rather complicated partitions of the domain of f into cells that are convex polytopes. When the architecture of the NN is fixed and the parameters are allowed to vary, the set of output functions of the NN is a parametrized nonlinear manifold. It is shown that this manifold has certain space-filling properties leading to an increased ability to approximate (better rate distortion) but at the expense of numerical stability. The space filling creates the challenge to the numerical method of finding best or good parameter choices when trying to approximate.

Contact phenomena arise in a variety of industrial process and engineering applications. For this reason, contact mechanics has attracted substantial attention from research communities. Mathematical problems from contact mechanics have been studied extensively for over half a century. Effort was initially focused on variational inequality formulations, and in the past ten years considerable effort has been devoted to contact problems in the form of hemivariational inequalities. This article surveys recent development in studies of hemivariational inequalities arising in contact mechanics. We focus on contact problems with elastic and viscoelastic materials, in the framework of linearized strain theory, with a particular emphasis on their numerical analysis. We begin by introducing three representative mathematical models which describe the contact between a deformable body in contact with a foundation, in static, history-dependent and dynamic cases. In weak formulations, the models we consider lead to various forms of hemivariational inequalities in which the unknown is either the displacement or the velocity field. Based on these examples, we introduce and study three abstract hemivariational inequalities for which we present existence and uniqueness results, together with convergence analysis and error estimates for numerical solutions. The results on the abstract hemivariational inequalities are general and can be applied to the study of a variety of problems in contact mechanics; in particular, they are applied to the three representative mathematical models. We present numerical simulation results giving numerical evidence on the theoretically predicted optimal convergence order; we also provide mechanical interpretations of simulation results.

In the last few decades the concepts of structure-preserving discretization, geometric integration and compatible discretizations have emerged as subfields in the numerical approximation of ordinary and partial differential equations. The article discusses certain selected topics within these areas; discretization techniques both in space and time are considered. Lie group integrators are discussed with particular focus on the application to partial differential equations, followed by a discussion of how time integrators can be designed to preserve first integrals in the differential equation using discrete gradients and discrete variational derivatives.

Lie group integrators depend crucially on fast and structure-preserving algorithms for computing matrix exponentials. Preservation of domain symmetries is of particular interest in the application of Lie group integrators to PDEs. The equivariance of linear operators and Fourier transforms on non-commutative groups is used to construct fast structure-preserving algorithms for computing exponentials. The theory of Weyl groups is employed in the construction of high-order spectral element discretizations, based on multivariate Chebyshev polynomials on triangles, simplexes and simplicial complexes.

The theory of mixed finite elements is developed in terms of special inverse systems of complexes of differential forms, where the inclusion of cells corresponds to pullback of forms. The theory covers, for instance, composite piecewise polynomial finite elements of variable order over polyhedral grids. Under natural algebraic and metric conditions, interpolators and smoothers are constructed, which commute with the exterior derivative and whose product is uniformly stable in Lebesgue spaces. As a consequence we obtain not only eigenpair approximation for the Hodge–Laplacian in mixed form, but also variants of Sobolev injections and translation estimates adapted to variational discretizations.