Past Numerical Analysis Group Internal Seminar

High-order tensor methods for solving both convex and nonconvex optimization problems have recently generated significant research interest, due in part to the natural way in which higher derivatives can be incorporated into adaptive regularization frameworks, leading to algorithms with optimal global rates of convergence and local rates that are faster than Newton's method. On each iteration, to find the next solution approximation, these methods require the unconstrained local minimization of a (potentially nonconvex) multivariate polynomial of degree higher than two, constructed using third-order (or higher) derivative information, and regularized by an appropriate power of the change in the iterates. Developing efficient techniques for the solution of such subproblems is currently, an ongoing topic of research,  and this talk addresses this question for the case of the third-order tensor subproblem.

In particular, we propose the CQR algorithmic framework, for minimizing a nonconvex Cubic multivariate polynomial with  Quartic Regularisation, by sequentially minimizing a sequence of local quadratic models that also incorporate both simple cubic and quartic terms. The role of the cubic term is to crudely approximate local tensor information, while the quartic one provides model regularization and controls progress. We provide necessary and sufficient optimality conditions that fully characterise the global minimizers of these cubic-quartic models. We then turn these conditions into secular equations that can be solved using nonlinear eigenvalue techniques. We show, using our optimality characterisations, that a CQR algorithmic variant has the optimal-order evaluation complexity of $O(\epsilon^{-3/2})$ when applied to minimizing our quartically-regularised cubic subproblem, which can be further improved in special cases.  We propose practical CQR variants that judiciously use local tensor information to construct the local cubic-quartic models. We test these variants numerically and observe them to be competitive with ARC and other subproblem solvers on typical instances and even superior on ill-conditioned subproblems with special structure.

Fourth-order elliptic problems arise in a variety of applications from thin plates to phase separation to liquid crystals. A conforming Galerkin discretization requires a finite dimensional subspace of $H^2$, which in turn means that conforming finite element subspaces are $C^1$-continuous. In contrast to standard $H^1$-conforming $C^0$ elements, $C^1$ elements, particularly those of high order, are less understood from a theoretical perspective and are not implemented in many existing finite element codes. In this talk, we address the implementation of the elements. In particular, we present algorithms that compute $C^1$ finite element approximations to fourth-order elliptic problems and which only require elements with at most $C^0$-continuity. We also discuss solvers for the resulting subproblems and illustrate the method on a number of representative test problems.

We present new high-order finite elements discretizing the $L^2$ de Rham complex on triangular and tetrahedral meshes. The finite elements discretize the same spaces as usual, but with different basis functions. They allow for fast linear solvers based on static condensation and space decomposition methods.

The new elements build upon the definition of degrees of freedom given by (Demkowicz et al., De Rham diagram for $hp$ finite element spaces. Comput.~Math.~Appl., 39(7-8):29--38, 2000.), and consist of integral moments on a symmetric reference simplex with respect to a numerically computed polynomial basis that is orthogonal in both the $L^2$- and $H(\mathrm{d})$-inner products ($\mathrm{d} \in \{\mathrm{grad}, \mathrm{curl}, \mathrm{div}\}$).

On the reference symmetric simplex, the resulting stiffness matrix has diagonal interior block, and does not couple together the interior and interface degrees of freedom. Thus, on the reference simplex, the Schur complement resulting from elimination of interior degrees of freedom is simply the interface block itself.

This sparsity is not preserved on arbitrary cells mapped from the reference cell. Nevertheless, the interior-interface coupling is weak because it is only induced by the geometric transformation. We devise a preconditioning strategy by neglecting the interior-interface coupling. We precondition the interface Schur complement with the interface block, and simply apply point-Jacobi to precondition the interior block.

The combination of this approach with a space decomposition method on small subdomains constructed around vertices, edges, and faces allows us to efficiently solve the canonical Riesz maps in $H^1$, $H(\mathrm{curl})$, and $H(\mathrm{div})$, at very high order. We empirically demonstrate iteration counts that are robust with respect to the polynomial degree.

Optimization problems constrained to a smooth manifold can be solved via the framework of Riemannian optimization. To that end, a geometrical description of the constraining manifold, e.g., tangent spaces, retractions, and cost function gradients, is required. In this talk, we present a novel approach that allows performing approximate Riemannian optimization based on a manifold learning technique, in cases where only a noiseless sample set of the cost function and the manifold’s intrinsic dimension are available.

Gaussian processes provide a powerful probabilistic kernel learning framework, which allows high-quality nonparametric learning via methods such as Gaussian process regression. Nevertheless, its learning phase requires unrealistic massive computations for large datasets. In this talk, we present a quadrature-based approach for scaling up Gaussian process regression via a low-rank approximation of the kernel matrix. The low-rank structure is utilized to achieve effective hyperparameter learning, training, and prediction. Our Gauss-Legendre features method is inspired by the well-known random Fourier features approach, which also builds low-rank approximations via numerical integration. However, our method is capable of generating high-quality kernel approximation using a number of features that is poly-logarithmic in the number of training points, while similar guarantees will require an amount that is at the very least linear in the number of training points when using random Fourier features. The utility of our method for learning with low-dimensional datasets is demonstrated using numerical experiments.

The proximal Galerkin finite element method is a high-order, nonlinear numerical method that preserves the geometric and algebraic structure of bound constraints in infinitedimensional function spaces. In this talk, we will introduce the proximal Galerkin method and apply it to solve free-boundary problems, enforce discrete maximum principles, and develop scalable, mesh-independent algorithms for optimal design. The proximal Galerkin framework is a natural consequence of the latent variable proximal point (LVPP) method, which is an stable and robust alternative to the interior point method that will also be introduced in this talk.

In particular, LVPP is a low-iteration complexity, infinite-dimensional optimization algorithm that may be viewed as having an adaptive barrier function that is updated with a new informative prior at each (outer loop) optimization iteration. One of the main benefits of this algorithm is witnessed when analyzing the classical obstacle problem. Therein, we find that the original variational inequality can be replaced by a sequence of semilinear partial differential equations (PDEs) that are readily discretized and solved with, e.g., high-order finite elements. Throughout the talk, we will arrive at several unexpected contributions that may be of independent interest. These include (1) a semilinear PDE we refer to as the entropic Poisson equation; (2) an algebraic/geometric connection between high-order positivity-preserving discretizations and an infinite-dimensional Lie group; and (3) a gradient-based, bound-preserving algorithm for two-field density-based topology optimization.

The complete latent variable proximal Galerkin methodology combines ideas from nonlinear programming, functional analysis, tropical algebra, and differential geometry and can potentially lead to new synergies among these areas as well as within variational and numerical analysis. This talk is based on [1].

Keywords: pointwise bound constraints, bound-preserving discretization, entropy regularization, proximal point

Mathematics Subject Classifications (2010): 49M37, 65K15, 65N30

References  [1] B. Keith, T.M. Surowiec. Proximal Galerkin: A structure-preserving finite element method for pointwise bound constraints arXiv preprint arXiv:2307.12444 2023.

Brown University Email address: @email

Simula Research Laboratory Email address: @email

Diffusion is one of the most common phenomenon in natural sciences and large part of applied mathematics have been interested in the tools to model it. Trying to study different types of diffusions, the mathematical ways to describe them and the numerical methods to simulate them is an appealing challenge, giving a wide range of applications. The aim of our work is the design of a finite-volume numerical scheme to model non-local diffusion given by the fractional Laplacian and to build numerical solutions for the Lévy-Fokker-Planck equation that involves it. Numerical methods for fractional diffusion have been indeed developed during the last few years and large part of the literature has been focused on finite element methods. Few results have been rather proposed for different techniques such as finite volumes.

 
We propose a new fractional Laplacian for bounded domains, which is expressed as a conservation law. This new approach is therefore particularly suitable for a finite volumes scheme and allows us also to prescribe no-flux boundary conditions explicitly. We enforce our new definition with a well-posedness theory for some cases to then capture with a good level of approximation the action of fractional Laplacian and its anomalous diffusion effect with our numerical scheme. The numerical solutions we get for the Lévy-Fokker-Planck equation resemble in fact the known analytical predictions and allow us to numerically explore properties of this equation and compute stationary states and long-time asymptotics.

The Ciarlet definition of a finite element has been used for many years to describe the requisite parts of a finite element. In that time, finite element theory and implementation have both developed and improved, which has left scope for a redefinition of the concept of a finite element. In this redefinition, we look to encapsulate some of the assumptions that have historically been required to complete Ciarlet’s definition, as well as incorporate more information, in particular relating to the symmetries of finite elements, using concepts from Group Theory. This talk will present the machinery of the proposed new definition, discuss its features and provide some examples of commonly used elements.

The PDE formulation of Mean Field Games (MFG) is described by nonlinear systems in which a Hamilton—Jacobi—Bellman (HJB) equation and a Kolmogorov—Fokker—Planck (KFP) equation are coupled. The advective term of the KFP equation involves a partial derivative of the Hamiltonian that is often assumed to be continuous. However, in many cases of practical interest, the underlying optimal control problem of the MFG may give rise to bang-bang controls, which typically lead to nondifferentiable Hamiltonians. In this talk we present results on the analysis and numerical approximation of second-order MFG systems for the general case of convex, Lipschitz, but possibly nondifferentiable Hamiltonians.
In particular, we propose a generalization of the MFG system as a Partial Differential Inclusion (PDI) based on interpreting the partial derivative of the Hamiltonian in terms of subdifferentials of convex functions.

We present theorems that guarantee the existence of unique weak solutions to MFG PDIs under a monotonicity condition similar to one that has been considered previously by Lasry & Lions. Moreover, we introduce a monotone finite element discretization of the weak formulation of MFG PDIs and prove the strong convergence of the approximations to the value function in the H¹-norm and the strong convergence of the approximations to the density function in L^q-norms. We conclude the talk with some numerical experiments involving non-smooth solutions. 

This is joint work with my supervisor Iain Smears. 

We introduce a sparse and very high order hp-finite element method for the weak form of the Helmholtz equation.  The domain may be a disk, an annulus, or a cylinder. The cells of the mesh are an innermost disk (omitted if the domain is an annulus) and concentric annuli.

We demonstrate the effectiveness of this method on PDEs with radial direction discontinuities in the coefficients and data. The discretization matrix is always symmetric and positive-definite in the positive-definite Helmholtz regime. Moreover, the Fourier modes decouple, reducing a two-dimensional PDE solve to a series of one-dimensional solves that may be computed in parallel, scaling with linear complexity. In the positive-definite case, we utilize the ADI method of Fortunato and Townsend to apply the method to a 3D cylinder with a quasi-optimal complexity solve.

The computation of matrix function is an important task appearing in many areas of scientific computing. We propose two algorithms, one for banded matrices and the other for general sparse matrices that approximate f(A) using matrix-vector products only. Our algorithms are inspired by the decay bound for the entries of f(A) when A is banded or sparse. We show its exponential convergence when A is banded or sufficiently sparse and we demonstrate its performance using numerical experiments.

For many stationary-state PDEs, solutions can be shown to satisfy certain key identities or structures, with physical interpretations such as the dissipation of energy. By reformulating these systems in terms of new auxiliary functions, finite-element models can ensure these structures also hold exactly for the numerical solutions. This approach is known to improve the solutions' accuracy and reliability.

In this talk, we extend this auxiliary function approach to the transient case through a finite-element-in-time interpretation. This allows us to develop novel structure-preserving timesteppers for various transient problems, including the Navier–Stokes and MHD equations, up to arbitrary order in time.

The Stokes–Onsager–Stefan–Maxwell (SOSM) equations model the flow of concentrated mixtures of distinct chemical species in a common thermodynamic phase. We derive a novel variational formulation of these nonlinear equations in which the species mass fluxes are treated as unknowns. This new formulation leads to a large class of high-order finite element schemes with desirable linear-algebraic properties. The schemes are provably convergent when applied to a linearization of the SOSM problem.

TBC

 Wide Neural Networks are well known for their Gaussian Process behaviour. Based upon this fact, an initialisation scheme for the weights and biases of a network preserving some geometrical properties of the input data is presented — The edge-of-chaos. This talk will introduce such a scheme before briefly mentioning a recent contribution related to the edge-of-chaos dynamics of wide randomly initialized low-rank feedforward networks. Formulae for the optimal weight and bias variances are extended from the full-rank to low-rank setting and are shown to follow from multiplicative scaling. The principle second order effect, the variance of the input-output Jacobian, is derived and shown to increase as the rank to width ratio decreases. These results inform practitioners how to randomly initialize feedforward networks with a reduced number of learnable parameters while in the same ambient dimension, allowing reductions in the computational cost and memory constraints of the associated network.

Placeholder entry; date+time TBC. 

Abstract for talk: In this talk, we present a systematic discretization of the Bernstein-Gelfand-Gelfand (BGG) diagrams and complexes over cubical meshes of arbitrary dimension via the use of tensor-product structures of one-dimensional piecewise-polynomial spaces, such as spline and finite element spaces. We demonstrate the construction of the Hessian, the elasticity, and div-div complexes as examples for our construction.

In this talk, we introduce Newton-MR variants for solving nonconvex optimization problems. Unlike the overwhelming majority of Newton-type methods, which rely on conjugate gradient method as the primary workhorse for their respective sub-problems, Newton-MR employs minimum residual (MINRES) method. With certain useful monotonicity properties of MINRES as well as its inherent ability to detect non-positive curvature directions as soon as they arise, we show that our algorithms come with desirable properties including the optimal first and second-order worst-case complexities. Numerical examples demonstrate the performance of our proposed algorithms.

We propose a nematic model for polyatomic gas, intending to study anisotropic phenomena. Such phenomena stem from the orientational degree of freedom associated with the rod-like molecules composing the gas. We adopt as a primer the Curitss-Boltzmann equation. The main difference with respect to Curtiss theory of hard convex body fluids is the fact that the model here presented takes into account the emergence of a nematic ordering. We will also derive from a kinetic point of view an energy functional similar to the Oseen-Frank energy. The application of the Noll-Coleman procedure to derive new expressions for the stress tensor and the couple-stress tensor will lead to a model capable of taking into account anisotropic effects caused by the emergence of a nematic ordering. In the near future, we hope to adopt finite-element discretisations together with multi-scale methods to simulate the integro-differential equation arising from our model.

Coupled differential equations generally present an important algebraic structure.
For example in the incompressible Navier-Stokes equations, the velocity is affected only by the selenoidal part of the applied force.
This structure can be translated naturally by the notion of complex.
One idea is then to exploit this complex structure at the discrete level in the creation of numerical methods.

The goal of the presentation is to expose the notion of complex by motivating its uses. 
We will present in more detail the creation of a scheme for the Navier-Stokes equations and study its properties.
 

The recent advancements in additive manufacturing have enabled the creation of lattice structures with intricate small-scale details. This has led to the need for new techniques in the field of topology optimization that can handle a vast number of design variables. Despite the efforts to develop multi-scale topology optimization techniques, their high computational cost has limited their application. To overcome this challenge, a new technique called dehomogenization has shown promising results in terms of performance and computational efficiency for optimizing compliance problems.

In this talk, we extend the application of the dehomogenization method to stress minimization problems, which are crucial in structural design. The method involves homogenizing the macroscopic response of a proposed family of microstructures. Next, the macroscopic structure is optimized using gradient-based methods while orienting the cells according to the principal stress components. The final step involves dehomogenization of the structure. The proposed methodology also considers singularities in the orientation field by incorporating singular functions in the dehomogenization process. The validity of the methodology is demonstrated through several numerical examples.

At the heart of all quasi-Newton methods is an update rule that enables us to gradually improve the Hessian approximation using the already available gradient evaluations. Theoretical results show that the global performance of optimization algorithms can be improved with higher-order derivatives. This motivates an investigation of generalizations of quasi-Newton update rules to obtain for example third derivatives (which are tensors) from Hessian evaluations. Our generalization is based on the observation that quasi-Newton updates are least-change updates satisfying the secant equation, with different methods using different norms to measure the size of the change. We present a full characterization for least-change updates in weighted Frobenius norms (satisfying an analogue of the secant equation) for derivatives of arbitrary order. Moreover, we establish convergence of the approximations to the true derivative under standard assumptions and explore the quality of the generated approximations in numerical experiments.

Sketch-and-precondition techniques are popular for solving large overdetermined least squares (LS) problems. This is when a right preconditioner is constructed from a randomized 'sketch' of the matrix. In this talk, we will see that the sketch-and-precondition technique is not numerically stable for ill-conditioned LS problems. We propose a modifciation: using an unpreconditioned iterative LS solver on the preconditioned matrix. Provided the condition number of A is smaller than the reciprocal of the unit round-off, we show that this modification ensures that the computed solution has a comparable backward error to the iterative LS solver applied to a well-conditioned matrix. Using smoothed analysis, we model floating-point rounding errors to provide a convincing argument that our modification is expected to compute a backward stable solution even for arbitrarily ill-conditioned LS problems.

We consider the problem of solving nonconvex quadratic optimization problems, potentially with additional integrality constraints on the variables.  Gurobi takes a branch-and-bound approach to solve such problems to global optimality, and in this talk we will review the three main algorithmic components that Gurobi uses:  Convex relaxations based on local linearization, globally valid cutting planes, and nonlinear local optimization heuristics.  We will explain how these parts play together, and discuss some of the implementation details.

The numerical solution of elliptic PDEs is often the most computationally intensive task in large-scale continuum mechanics simulations.  High-order finite element methods can efficiently exploit modern parallel hardware while offering very rapid convergence properties.  As the polynomial degree is increased, the efficient solution of such PDEs becomes difficult. In this talk we introduce preconditioners for high-order discretizations. We build upon the pioneering work of Pavarino, who proved in 1993 that the additive Schwarz method with vertex patches and a low-order coarse space gives a  solver for symmetric and coercive problems that is robust to the polynomial degree. 

However, for very high polynomial degrees it is not feasible to assemble or factorize the matrices for each vertex patch, as the patch matrices contain dense blocks, which couple together all degrees of freedom within a cell. The central novelty of the preconditioners we develop is that they have the same time and space complexity as sum-factorized operator application on unstructured meshes of tensor-product cells, i.e., we can solve $A x=b$ with the same complexity as evaluating $b-A x$. Our solver relies on new finite elements for the de Rham complex that enable the blocks in the stiffness matrix corresponding to the cell interiors to become diagonal for scalar PDEs or block diagonal for vector-valued PDEs.  With these new elements, the patch problems are as sparse as a low-order finite difference discretization, while having a sparser Cholesky factorization. In the non-separable case, themethod can be applied as a preconditioner by approximating the problem with a separable surrogate.  Through the careful use of incomplete factorizations and choice of space decomposition we achieve optimal fill-in in the patch factors, ultimately allowing for optimal-complexity storage and computational cost across the setup and solution stages.

We demonstrate the approach by solving the Riesz maps of $H^1$, $H(\mathrm{curl})$, and $H(\mathrm{div})$ in three dimensions at $p = 15$.

We investigate the theoretical properties of subsampling and hashing as tools for approximate Euclidean norm-preserving embeddings for vectors with (unknown) additive Gaussian noises. Such embeddings are called Johnson-Lindenstrauss embeddings due to their celebrated lemma. Previous work shows that as sparse embeddings, if a comparable embedding dimension to the Gaussian matrices is required, the success of subsampling and hashing closely depends on the $l_\infty$ to $l_2$ ratios of the vectors to be mapped. This paper shows that the presence of noise removes such constrain in high-dimensions; in other words, sparse embeddings such as subsampling and hashing with comparable embedding dimensions to dense embeddings have similar norm-preserving dimensionality-reduction properties, regardless of the $l_\infty$ to $l_2$ ratios of the vectors to be mapped. The key idea in our result is that the noise should be treated as information to be exploited, not simply a nuisance to be removed. Numerical illustrations show better performances of sparse embeddings in the presence of noise.

In this talk we are presenting a new approach for compatible finite element discretisations for atmospheric flows on a terrain following mesh. In classical compatible finite element discretisations, the H(div)-velocity space involves the application of Piola transforms when mapping from a reference element to the physical element in order to guarantee normal continuity. In the case of a terrain following mesh, this causes an undesired coupling of the horizontal and vertical velocity components. We are proposing a new finite element space, that drops the Piola transform. For solving the equations we introduce a hybridisable formulation with trace variables supported on horizontal cell faces in order to enforce the normal continuity of the velocity in the solution. Alongside the discrete formulation for various fluid equations we discuss solver approaches that are compatible with them and present our latest numerical results.

While scattering problems are posed on unbounded domains, volumetric discretizations typically require truncating the domain at a finite distance, closing the system with some sort of boundary condition.  These conditions typically suffer from some deficiency, such as perturbing the boundary value problem to be solved or changing the character of the operator so that the discrete system is difficult to solve with iterative methods.

We introduce a new technique for the Helmholtz problem, based on using the Green formula representation of the solution at the artificial boundary.  Finite element discretization of the resulting system gives optimal convergence estimates.  The resulting algebraic system can be solved effectively with a matrix-free GMRES implementation, preconditioned with the local part of the operator.  Extensions to the Morse-Ingard problem, a coupled system of pressure/temperature equations arising in modeling trace gas sensors, will also be given.

In this talk I will present a new class of algorithms for separable nonlinear inverse problems based on inexact Krylov methods. In particular, I will focus in semi-blind deblurring applications. In this setting, inexactness stems from the uncertainty in the parameters defining the blur, which are computed throughout the iterations. After giving a brief overview of the theoretical properties of these methods, as well as strategies to monitor the amount of inexactness that can be tolerated, the performance of the algorithms will be shown through numerical examples. This is joint work with Silvia Gazzola (University of Bath).

Rational functions are able to approximate functions containing branch point singularities with a root-exponential convergence rate. These appear for example in the solution of boundary value problems on domains containing corners or edges. Results from Newman in 1964 indicate that the poles of the optimal rational approximant are exponentially clustered near the branch point singularities. Trefethen and collaborators use this knowledge to linearize the approximation problem by fixing the poles in advance, giving rise to the Lightning approximation. The resulting approximation set is however highly ill-conditioned, which raises the question of stability. We show that augmenting the approximation set with polynomial terms greatly improves stability. This observation leads to a  decoupling of the approximation problem into two regimes, related to the singular and the smooth behaviour of the function. In addition, adding polynomial terms to the approximation set can result in a significant increase in convergence speed. The convergence rate is however very sensitive to the speed at which the clustered poles approach the singularity.

Most algorithms for computing a matrix function $f(A)$ are based on finding a rational (or polynomial) approximant $r(A)≈f(A)$ to the scalar function on the spectrum of $A$. These functions are often in a composite form, that is, $f(z)≈r(z)=r_k(⋯r_2(r_1(z)))$ (where $k$ is the number of compositions, which is often the iteration count, and proportional to the computational cost); this way $r$ is a rational function whose degree grows exponentially in $k$. I will review algorithms that fall into this category and highlight the remarkable power of composite (rational) functions.

An occupational hazard of mathematicians is the investigation of objects that are "optimal" in a mathematically precise sense, yet may be far from optimal in practice. This talk will discuss an extreme example of this effect: Gauss-Hermite quadrature on the real line. For large numbers of quadrature points, Gauss-Hermite quadrature is a very poor method of integration, much less efficient than simply truncating the interval and applying Gauss-Legendre quadrature or the periodic trapezoidal rule. We will present a theorem quantifying this difference and explain where the standard notion of optimality has failed.

Motivated by the physical relevance of many Hodge Laplace-type PDEs from the finite element exterior calculus, we analyse the Hodge Laplacian boundary value problem arising from the strain space in the linear elasticity complex, an exact sequence of function spaces naturally arising in several areas of continuum mechanics. We propose a discretisation based on the adaptation of discontinuous Galerkin FEM for the incompatibility operator $\mathrm{inc} := \mathrm{rot}\circ\mathrm{rot}$, using the symmetric-tensor-valued Regge finite element to discretise  the strain field; via the 'Regge calculus', this element has already been successfully applied to discretise another metric tensor, namely that arising in general relativity. Of central interest is the characterisation of the associated Sobolev space $H(\mathrm{inc};\mathbb{R}^{d\times d}_{\mathrm{sym}})$. Building on the pioneering work of van Goethem and coauthors, we also discuss promising connections between functional analysis of the $\mathrm{inc}$ operator and Kröner's theory of intrinsic elasticity in the presence of defects.

This is based on ongoing work with Dr Kaibo Hu.

Regularization of linear least squares problems is necessary in a variety of contexts. However, the optimal regularization parameter is usually unknown a priori and is often to be determined in an ad hoc manner, which may involve solving the problem for multiple regularization parameters. In this talk, we will discuss three randomized algorithms, building on the sketch-and-precondition framework in randomized numerical linear algebra (RNLA), to efficiently solve this set of problems. In particular, we consider preconditioners for a set of Tikhonov regularization problems to be solved iteratively. The first algorithm is a Cholesky-based algorithm employing a single sketch for multiple parameters; the second algorithm is SVD-based and improves the computational complexity by requiring a single decomposition of the sketch for multiple parameters. Finally, we introduce an algorithm capable of exploiting low-rank structure (specifically, low statistical dimension), requiring a single sketch and a single decomposition to compute multiple preconditioners with low-rank structure. This algorithm avoids the Gram matrix, resulting in improved stability as compared to related work.

Time-optimal control can be used to improve driving efficiency for autonomous
vehicles and it enables us explore vehicle and driver behaviour in extreme
situations. Due to the computational cost and limited scope of classical
optimal control methods we have seen new interest in applying reinforcement
learning algorithms to autonomous driving tasks.
In this talk we present methods for translating time-optimal vehicle control
problems into reinforcement learning environments. For this translation we
construct a sequence of environments, starting from the closest representation
of our optimisation problem, and gradually improve the environments reward
signal and feature quality. The trained agents we obtain are able to generalise
across different race tracks and obtain near optimal solutions, which can then
be used to speed up the solution of classical time-optimal control problems.

Many problems in engineering can be understood as controlling the bifurcation structure of a given device. For example, one may wish to delay the onset of instability, or bring forward a bifurcation to enable rapid switching between states. In this talk, we will describe a numerical technique for controlling the bifurcation diagram of a nonlinear partial differential equation by varying the shape of the domain or a parameter in the equation. Our aim is to delay or advance a given branch point to a target parameter value. The algorithm consists of solving an optimization problem constrained by an augmented system of equations that characterize the location of the branch points. The flexibility and robustness of the method also allow us to advance or delay a Hopf bifurcation to a target value of the bifurcation parameter, as well as controlling the oscillation frequency. We will apply this technique on systems arising from biology, fluid dynamics, and engineering, such as the FitzHugh-Nagumo model, Navier-Stokes, and hyperelasticity equations.

This talk is about the mathematics behind an artistic project focusing on the vibrations of soap films.

We investigate the maximum distance of a rank-two tensor to rank-one tensors. An equivalent problem is given by the minimal ratio of spectral and Frobenius norm of a tensor. For matrices the distance of a rank k matrix to a rank r matrices is determined by its singular values, but since there is a lack of a fitting analog of the singular value decomposition for tensors, this question is more difficult in the regime of tensors.
 

iii) identify a special class of correlated compressors based on the idea of random permutations, for which we coin the term PermK. The use of this technique results in the strict improvement on the previous MARINA rate. In the low Hessian variance regime, the improvement can be as large as √n, when d > n, and 1 + √d/n, when n<=d, where n is the number of workers and d is the number of parameters describing the model we are learning.

We introduce meta-factorization, a theory that describes matrix decompositions as solutions of linear matrix equations: the projector and the reconstruction equation. Meta-factorization reconstructs known factorizations, reveals their internal structures, and allows for introducing modifications, as illustrated with SVD, QR, and UTV factorizations. The prospect of meta-factorization also provides insights into computational aspects of generalized matrix inverses and randomized linear algebra algorithms. The relations between the Moore-Penrose pseudoinverse, generalized Nyström method, and the CUR decomposition are revealed here as an illustration. Finally, meta-factorization offers hints on the structure of new factorizations and provides the potential of creating them. 

Mass transfer in multicomponent systems occurs through convection and diffusion. For a viscous Newtonian flow, convection may be modelled using the Navier–Stokes equations, whereas the diffusion of multiple species within a common phase may be described by the generalised Onsager–Stefan–Maxwell equations. In this talk we present a novel finite element formulation which fully couples convection and diffusion with these equations. In the regime of vanishing Reynolds number, we use the principles of linear irreversible dynamics to formulate a saddle point system which leads to a stable formulation and a convergent discretisation. The wide scope of applications for this novel numerical method is illustrated by considering transport of oxygen through the lungs, gas separation processes, mixing of water and methanol and salt transport in electrolytes.

The placenta is the essential organ of maternal-fetal interactions, where nutrient, oxygen, and waste exchange occur. In recent studies, differences in the morphology of the placental chorionic surface vascular network (PCSVN) have been associated with developmental disorders such as autism. This suggests that the PCSVN could potentially serve as a biomarker for the early diagnosis and treatment of autism. Studying PCSVN features in large cohorts requires a reliable and automated mechanism to extract the vascular networks. In this talk, we present a method for PCSVN extraction. Our algorithm builds upon a directional multiscale mathematical framework based on a combination of shearlets and Laplacian eigenmaps and can isolate vessels with high success in high-contrast images such as those produced in CT scans. 

TBA.

How can one integrate singular functions over fractals? And why would one want to do this? In this talk I will present a general approach to numerical quadrature on the compact attractor of an iterated function system of contracting similarities, where integration is with respect to the relevant Hausdorff measure. For certain singular integrands of logarithmic or algebraic type the self-similarity of the integration domain can be exploited to express the singular integral exactly in terms of regular integrals that can be approximated using standard techniques. As an application we show how this approach, combined with a singularity-subtraction technique, can be used to accurately evaluate the singular double integrals that arise in Hausdorff-measure Galerkin boundary element methods for acoustic wave scattering by fractal screens. This is joint work with Andrew Gibbs (UCL) and Andrea Moiola (Pavia).

Many applications involve non-Euclidean data, where exploiting Riemannian geometry can deliver algorithms that are computationally superior to standard nonlinear programming approaches. This observation has resulted in an increasing interest in Riemannian methods in the optimization and machine learning community. In this talk, we consider the problem of optimizing a function on a Riemannian manifold subject to convex constraints. We introduce Riemannian Frank-Wolfe (RFW) methods, a class of projection-free algorithms for constrained geodesically convex optimization. To understand the algorithm’s efficiency, we discuss (1) its iteration complexity, (2) the complexity of computing the Riemannian gradient and (3) the complexity of the Riemannian “linear” oracle (RLO), a crucial subroutine at the heart of the algorithm. We complement our theoretical results with an empirical comparison of RFW against state-of-the-art Riemannian optimization methods. Joint work with Suvrit Sra (MIT).