Numerical Analysis and Scientific Computing Seminars 2011/12
Semester Two (Spring 2012)

03 Feb
2012
New generation of numerical algorithms via tensortrain representations of data
Eugene Tyrtyshnikov (Russian Academy of Sciences, Moscow)
3.00 pm  Frank Adams Room 1, Alan Turing BuildingAbstract (click to view)In 20th century tensors were used chiefly as a language or a descriptive means. In 21st century they get to play a new role as a base for numerical algorithms for the data in really many (as well as in few) dimensions.
The key is a new tensor decomposition called the tensortrain (TT) decomposition and basic TT algorithms appeared only in 2009: TT rounding, TT cross interpolation, and TT wavelet transorms. We discuss these new findings and the related perspectives for the development of numerical analysis. 
10 Feb
2012
Identifying Hopf Bifurcations in Models of Incompressible Flow
Alastair Spence (University of Bath)
3.00 pm  Frank Adams Room 1, Alan Turing BuildingAbstract (click to view)The identification of instability in largescale dynamical systems caused by Hopf bifurcation is difficult because of the problem of computing the rightmost pair of complex eigenvalues of large sparse generalised eigenvalue problems. This talk describes a method developed in [Meerbergen & Spence, SIMAX (2010), pp. 19821999] that avoids this computation, instead performing an inverse iteration for a certain set of real eigenvalues that requires the solution of a largescale Lyapunov equation at each iteration. The method is tested on two challenging problems arising from fluid dynamics.

17 Feb
2012
Condition Numbers in Optimisation and Sparse Recovery
Martin Lotz (University of Edinburgh)
3.00 pm  Frank Adams Room 1, Alan Turing BuildingAbstract (click to view)The running time of many iterative algorithms in numerical analysis and optimisation can be analysed in terms of condition numbers. Moreover, the condition numbers are often intimately related to geometric features of the problem. This talk surveys recent results on the probabilistic analysis of condition numbers. It is also shown how the analysis of condition numbers for convex optimisation problems leads to new approaches to fundamental questions related to the active field of compressed sensing.

24 Feb
2012
Software abstractions for manycore software engineering
Paul Kelly (Imperial College)
3.00 pm  Frank Adams Room 1, Alan Turing BuildingAbstract (click to view)What is the right code to generate, for a given hardware platform? How does this change as problem parameters change? This talk presents some recent workinprogress exploring domainspecific languages and active libraries as a way to automate code generation for multicore and manycore platforms, and to capture the space of alternative implementation choices at a higher level than a compiler for a generalpurpose language can. I will illustrate the potential for this idea by looking at our recent work on unstructuredmesh fluid dynamics, and finite element methods. By choosing the abstraction carefully, we can capture design choices far beyond what a conventional compiler can do  and, in the extreme, engage the users in selecting algorithms and numerical methods that match the capabilities of the underlying hardware to meet endtoend objectives for solution quality.
This talk presents joint work with David Ham, Spencer Sherwin and Chris Cantwell at Imperial, and Mike Giles, and Gihan Mudalige at Oxford. 
16 Mar
2012
Computable bounds for the error in finite element approximations
Richard Rankin (University of Strathclyde)
3.00 pm  Frank Adams Room 1, Alan Turing BuildingAbstract (click to view)Finite element approximation schemes can be used to obtain approximate solutions to partial differential equations. When using a numerical method to approximate the solution to any problem it is important to be able to say how accurate the approximation obtained is. I will discuss how computable estimators which provide guaranteed upper bounds on the energy norm of the error in finite element approximations can be obtained. I will also discuss the parallel implementation of these estimators.

23 Mar
2012
Solution of illposed inverse problems pertaining to signal restoration
Rosemary Renaut (Arizona State University, USA)
3.00 pm  Frank Adams Room 1, Alan Turing BuildingAbstract (click to view)In this talk I review the use of the spectral decomposition for understanding the solution of illposed inverse problems. It is immediate to see that regularization is needed in order to find stable solutions. These solutions, however, do not typically allow reconstruction of signal features such as edges. Generalized regularization assists but is still insufficient and methods of total variation are commonly suggested as an alternative. In the talk I consider application of standard approaches from Tikhonov regularization for finding appropriate regularization parameters in the total variation augmented Lagrangian implementations. The analysis of the solution uses the SVD and GSVD, and assists with understanding how to pick relevant regularization parameters for the total variation implementation. Areas for future research will be considered.

27 Apr
2012
Exponential integrators and porous media flow
Gabriel Lord (HeriotWatt University)
3.00 pm  Frank Adams Room 1, Alan Turing BuildingAbstract (click to view)We consider exponential based integrators for porous media flow in two and three dimensions. These methods are based on solving the linear problem exactly and exploit efficient methods for computing large matrix exponentials. For Krylov based methods we briefly discuss an implementation for timestepping that reuses the Arnoldi iterations, extending ideas in [1].
We are particularly interested in the effects of time dependent stochastic forcing. The forcing is taken to model unknown movement between trapped and flowing solutes, unknown small scale variations in the reaction term or time dependent changes in permeability. We examine convergence and efficiency of new schemes for the stochastic PDEs for both additive and multiplicative noise and present numerical results for realistic flow type problems.
References
[1] E. J. Carr, T. J. Moroney and I. W. Turner, Efficient simulation of unsaturated flow using exponential time integration, Applied Math. and Computation, 215 (2011), 65876596. 
25 May
2012
First Order kth Moment Finite Element Analysis
of Nonlinear Operator Equations with Stochastic Data
Alexey Chernov (Hausdorff Center for Mathematics, Bonn, Germany)
3.00  Frank Adams Room 1, Alan Turing BuildingAbstract (click to view)We develop and analyze a class of efficient algorithms for uncertainty quantification of nonlinear operator equations. The algorithms are based on sparse Galerkin discretizations of tensorized linearizations at nominal parameters. Specifically, we analyze a class of abstract nonlinear, parametric operator equations
J(α,u)=0
for random parameters α with realizations in a neighborhood of a nominal parameter α_{0}. Under some structural assumptions on the parameter dependence, random parameters α(ω) = α_{0} + r(ω) are shown to imply a unique random solution u(ω)=S(α(ω)).
We derive an operator equation for the deterministic computation of kth order statistical moments of the solution fluctuations u(ω)  S(α_{0}), provided that statistical moments of the random parameter perturbation r(ω) are known. This formulation is based on the linearization in a neighborhood of a nominal parameter α_{0}. We provide precise bounds for the linearization error by means of a NewtonKantorovichtype theorem.
We present a sparse tensor Galerkin discretization for the tensorized first order perturbation equation and deduce that sparse tensor Galerkin discretizations of this equation converge in accuracy vs. complexity which equals, up to logarithmic terms, that of the Galerkin discretization of a single instance of the mean field problem.
Joint work with Christoph Schwab (Seminar for Applied Mathematics, ETH Zürich)
References
[1] A. Chernov, Abstract sensitivity analysis for nonlinear equations and applications, In Karl Kunisch, Günther Of, and Olaf Steinbach, editors, Numerical Mathematics and Advanced Applications, Proceedings of ENUMATH 2007, Graz, Austria, pages 407414. Springer, September 2008.
[2] A. Chernov and C. Schwab, First order kth moment Finite Element analysis of nonlinear operator equations with stochastic data, Mathematics of Computation, accepted, also available as Preprint 2011b06, Hausdorff Research Institute for Mathematics, University of Bonn, August 2011, http://131.220.77.52/files/preprints/AnalysisandNumerics/ChernovSchwab_FOkM.pdf 
31 May
2012
Stabilization for conservation laws
Alexandre Ern (Université ParisEst, France)
3.00 pm  Frank Adams Room 1, Alan Turing BuildingAbstract (click to view)A modification of the edge stabilization technique is proposed to improve the behavior of the method when solving conservation equations with nonsmooth data and/or nonsmooth solutions. The key ingredient is tempering the edge stabilization in regions of large gradients through appropriate weights. The new method is shown to preserve the convergence properties of the original method on smooth solutions and numerical tests indicate that it performs better on nonsmooth solutions.
Semester One (Autumn 2011)

21 Sept
2011
Fast Algorithms for Approximating the Distance to
Instability and Hinfinity Norm
Mert Gurbuzbalaban (Courant Institute, New York, USA)
2.00  Frank Adams Room 1, Alan Turing Building
(Please note unusual day and time!)Abstract (click to view)The stability radius of an n x n square matrix A (or distance to instability) is a wellknown measure of the stability of the linear dynamical system dx/dt = Ax. Existing techniques compute this quantity accurately but the cost is multiple SVDs of order n, which makes the method suitable to middle size problems. We present a novel approach based on a Newton iteration applied to the pseudospectral abscissa, whose implementation is obtained by a repeated computation of the spectral abscissa of a sequence of matrices. Such an approach turns out to be particularly well suited for large sparse matrices. We will also discuss the extensions of this approach to the more general problem of computing the Hinfinity norm of transfer matrix.

29 Sept
2011
Weak Mass Conservation in Discretizations of the Incompressible
NavierStokes Equations
Alexander Linke (WIAS, Berlin, Germany)
3.15  Frank Adams Room 1, Alan Turing Building
(Please note unusual day and time!)Abstract (click to view)In this presentation, the instability of weak mass conservation in finite element and finite volume discretizations of the incompressible NavierStokes equations will be investigated. The lack of discrete L^{2}orthogonality of divergencefree and rotationfree functions will appear as a common problem of many methods which can trigger spurious numerical oscillations in the discrete velocities. Some examples from applications will be presented where indeed weak mass conservation is the most important problem in the numerical computation. Further, different approaches of mitigating or removing weak mass conservation will be presented, and the relation of these approaches to each other will be analyzed.

28 Oct
2011
Some Recent Advances in Computational Electromagnetics:
Solutions of LargeScale Problems Discretized With Hundreds of Millions of Unknowns
Ozgur Ergul (University of Strathclyde)
3.00  Frank Adams Room 1, Alan Turing BuildingAbstract (click to view)Recently, there have been many efforts for more and more realistic simulations of electromagnetic phenomena, which can be found across the entire spectrum of healthcare, wireless communications, aerospace industry, and energy. Computational simulations of those practical scenarios can provide the essential information on electronic devices and new designs before their actual realizations, preventing the waste of sources and time involved in building prototypes and carrying out laboratory tests. However, accurate simulations of reallife problems are challenging and require a welldesigned combination of different components from diverse areas, such as wave theory, computer modeling, integral equations, discretization, numerical methods, iterative techniques, preconditioning, fast algorithms, parallelization, and highperformance computing. Specifically, in electromagnetics, accurate modeling of realistic structures often leads to dense matrix equations involving millions of unknowns, which cannot be solved easily, even when using the most powerful computers in today's technology. Therefore, it becomes essential to develop special acceleration algorithms and utilize advanced parallelization techniques in order to solve these realistic problems efficiently and accurately. This talk will be on the development and implementation of a parallel solver and its application to largescale problems discretized with hundreds of millions of unknowns. Using a fullwave algorithm, namely, the multilevel fast multipole algorithm, and its parallelization via the hierarchical partitioning strategy, such extremely large problems can be solved rigorously on moderate computers with distributedmemory architectures.

11 Nov
2011
Efficient hpbem on curved surfaces with numerical quadrature
Matthias Maischak (Brunel University)
3.00  Frank Adams Room 1, Alan Turing BuildingAbstract (click to view)The hpmethod on geometrical refined meshes for 3DBEM leads to exponential fast convergence on polyhedral domains [4,5] with respect to the number of unknowns, as well as the pversion for 3dBEM on smooth domains. Using numerical quadrature, exponential fast convergence can be still achieved for curved surfaces [2,6], if the surface parametrisation is relatively cheap. Here we will make use of a general surface parametrisation, cf. [3], and we will allow nonmatching grids, cf. [1], to simplify the connection of separately meshed parts of the surface. Emphasis is given to balancing the computational costs for evaluation of the surface parametrisation, incorporating the constraints on the spline space due to hanging nodes, with the exponential fast convergence of the method.
We will present several examples where we obtain exponential fast convergence on complicated geometries.
References
[1] L. Demkowicz, K. Gerdes, C. Schwab, A. Bajer, and T. Walsh, HP90: A general and flexible Fortran 90 hpFE code, Comput. Visual. Sci., 1 (1998), pp. 145163.
[2] S. Erichsen and S. A. Sauter, Efficient automatic quadrature in 3d Galerkin BEM, Computer Methods in Applied Mechanics and Engineering, 157 (1998), pp. 215224.
[3] H. Harbrecht and R. Schneider, Wavelet Galerkin schemes for boundary integral equations, SIAM J. Sci. Comput., 27 (2006), pp. 13471370.
[4] N. Heuer, M. Maischak, and E. P. Stephan, Exponential convergence of the hpversion for the boundary element method on open surfaces, Numer. Math., 83 (1999), pp. 641666.
[5] H. Holm, M. Maischak, and E. P. Stephan, Exponential convergence of the hp version bem for mixed bvp's on polyhedrons, Math. Meth. Appl. Sci., 31 (2008), pp. 20692093.
[6] S. A. Sauter and C. Schwab, Quadrature for hpGalerkin BEM in R3, Numer. Math., 78 (1997), pp. 211258. 
18 Nov
2011
Mathematical techniques for Computer Architecture Design
George A. Constantinides (Imperial College)
3.00  Frank Adams Room 1, Alan Turing BuildingAbstract (click to view)This talk will introduce some of our work in the Circuits and Systems group at Imperial College, with emphasis on research spanning the broad spectrum of Numerical Analysis, HighPerformance Computing, Computer Architecture and Electronics Design. After a brief overview of the current directions being pursued in computer architecture, we will look at some of the challenges in mapping algorithms to highly flexible parallel architectures. The tension between precision and performance will be investigated from a hardware perspective, and placed in an abstract context to highlight our current approaches to addressing this question, based on real algebraic geometry. We will also look at the use of parametric optimisation for designing efficient memory systems for highperformance architectures. Finally, we will discuss some open questions in the area.

25 Nov
2011
Shape Analysis
Stephen Marsland (Massey University, New Zealand)
3.00  Frank Adams Room 1, Alan Turing BuildingAbstract (click to view)I will give an overview of the field of shape analysis, which looks for the mathematical, computational, and statistical tools required to recognise and interpret shapes in images. It has many uses, such as in medical image analysis, where identifying shape and shape change can assist in the diagnosis of disease. I will describe methods of shape analysis from the viewpoint of generalised Euler equations and geometric mechanics.

02 Dec
2011
On the numerical computation of the modified energy
Jitse Niesen (University of Leeds)
3.00  Frank Adams Room 1, Alan Turing BuildingAbstract (click to view)The numerical solution of an ordinary differential equation is the exact solution of a slightly different equation, called the modified equation. If the original differential equation is Hamiltonian and the numerical method is symplectic, then the modified equation is also Hamiltonian; the corresponding energy is called the modified energy.
The modified energy can be computed analytically but this is a complicated procedure. Skeel and coworkers found a remarkable method for computing the modified energy numerically and showed that the results shed light on the behaviour of the numerical integrator. We revisit their method and simplify it by using Richardson interpolation allowing us to attain very high order. This allows us to capture the exponentially small terms in backward error analysis causing the energy to drift.
This is joint work with Per Christian Moan.
Further information
For further information please contact the seminar organiser Alex Bespalov.