Morning Starter Coffee

wherein there will be coffee and tea get started

Welcome

A quick welcome to go over useful information about the workshop.

Monolithic multigrid solvers for high-order discretizations of saddle-point systems

Recent years have seen substantial interest in the development of high-order discretizations for many saddle-point systems, including the Stokes and Navier-Stokes equations and the equations of magnetohydrodynamics (MHD). While monolithic multigrid is well-known to be an h-robust preconditioner for corresponding low-order discretizations, its extension to higher-order discretizations has not been well-studied. In this talk, we discuss how to successfully extend these approaches from low-order to high-order discretizations.

Parallel preconditioning for time-harmonic waves from Toeplitz structured sweeping

We consider the potential for parallelism within an approach to precondition problems involving high-frequency time-harmonic waves. We take inspiration from domain decomposition strategies known as sweeping methods, which have gained notable interest for their ability to yield nearly-linear asymptotic complexity and which can also be favourable for high-frequency problems. While successful approaches exist, such as those based on higher-order interface conditions, perfectly matched layers (PMLs), or complex tracking of wave fronts, they can often be quite involved or tricky to implement and inherently sequential in nature. We introduce here parallelisable versions of simple sweeping techniques for the Helmholtz equation that allow for each subdomain to be solved in parallel with minimal impact on the overall convergence.

Coffee Break

wherein there will be coffee, tea, and small bites to eat

A quantum solver for the Navier-Stokes equations

We present an end-to-end quantum algorithm for simulating nonlinear dissipative PDEs. Our algorithm approximates the expected value of any correlation function that depends on $O(1)$ variables, with rigorous bounds on the approximation error. The runtime scales polylogarithmically in the number of degrees of freedom, and linearly with the evolution time and nonlinearity strength. The considered simulation problem is BQP-complete, providing strong evidence for a quantum advantage in the simulation of turbulent fluid dynamics.

Block Alpha-Circulant Preconditioners for All-at-Once Diffusion-Based Covariance Operators

Diffusion-based correlation operators are commonly-used in ocean data assimilation to weight the contribution of prior information. Applying these operators can require a large number of serial matrix-vector products. A parallel-in-time formulation removes this requirement, and offers the opportunity to exploit modern computer architectures. High quality preconditioners for the parallel-in-time approach are well-known, but impossible to apply in practice for the high-dimensional problems that occur in oceanography. In this talk we consider a nested preconditioning approach which retains many of the beneficial properties of the ideal analytic preconditioner while remaining affordable in terms of memory and computational resource.

Lunch

Lunch catered from Sweet Cecily

Preconditioning Random Feature Methods and Hybrid Iterative Methods

This talk examines two complementary sources of failure in machine-learning-based PDE solvers: nonphysical stagnation in hybrid iterative schemes and ill-conditioning in domain-decomposed random feature discretizations. For deep-learning-based hybrid iterative methods, I will show that reliability depends crucially on the interaction between training objectives and solver updates, and that physics-aware Anderson acceleration can recover convergence by minimizing the physical residual. For random feature methods, I will present a local RRQR filtering and right-preconditioning framework that removes redundant basis functions, preserves sparsity, and yields sharp conditioning improvements together with major LSQR speedups. The overarching message is that robust AI-PDE solvers require not only expressive models, but also numerically consistent iteration and linear algebra.

Preconditioning for Singularly Perturbed PDEs

We'll consider the solution of linear systems arising from the numerical solution of singularly perturbed partial differential equations. Such problems feature a perturbation parameter, and the literature is awash with studies of so-called parameter robust discretizations, whose accuracy is proven to be independent of the perturbation parameter. That is, for a fixed number of degrees of freedom, there is a guaranteed error bound, no matter how small the perturbation parameter. However, for such methods to be truly robust, the cost of solving the resulting linear system should also be independent of the perturbation parameter. I'll give a overview of the challenges faced by both direct and iterative solvers for such systems. Moreover, I'll discuss how, with careful preconditioning, parameter robust solvers can be designed.

Coffee Break

wherein there will be coffee, tea, and small bites to eat

Self-adjointness and complex eigenvalues of real matrices

Motivated by trying to understand the convergence of the GMRES Krylov subspace iterative method for the solution of linear systems of equations with a real nonsymmetric coefficient matrix, we present our recent theoretical results relating to complex eigenvalues. Specifically, because every real matrix is self-adjoint in a symmetric bilinear form (indefinite inner product), we generalise the well known result that a real symmetric matrix (in fact any real matrix that is self-adjoint in a definite inner product) has all real eigenvalues, by proving that the number of complex conjugate eigenpairs is bounded by the inertia of any real symmetric matrix defining such an indefinite inner product in which the matrix is self-adjoint. If such a matrix has few negative eigenvalues, then there are necessarily few complex eigenpairs. Simple numerical examples will be given which highlight these elementary (but seemingly new) mathematical results.

Freedom to Prescribe Convergence Behavior of GMRES on Tridiagonal Toeplitz Systems

Discretizing PDEs leads to linear systems with large, sparse coefficient matrices. When linear, constant-coefficient PDEs with Dirichlet boundary conditions are discretized on uniform meshes, one can obtain Toeplitz, multilevel Toeplitz and/or block Toeplitz systems. Toeplitz matrices have constant diagonals, and multilevel and block Toeplitz matrices have related structures, that can be exploited to speed up GMRES, and aid convergence analysis. Such systems are widely solved by Krylov subspace methods. When A is nonsymmetric, however, GMRES convergence behavior is complicated to describe. It is proven through APS parametrization of A that any nonoincreasing convergence curve is possible for GMRES by constructing a linear system of prescribed nonzero eigenvalues with a given convergence curve. Nevertheless, a reconstructed matrix A in this case does not preserve the Toeplitz structure in general. We show, also through APS parametrization of A, that it is generally hard to prescribe convergence behavior of GMRES on tridiagonal Toeplitz systems.

Collaboration Time

Morning Starter Coffee

wherein there will be coffee and tea get started

Dynamic updating and learning for preconditioned iterative solvers

Preconditioners are essential for the fast convergence of iterative linear solvers, but the computation of a good preconditioner can be expensive (for example, algebraic multigrid preconditioners). For matrix-free solvers, computing a good preconditioner is even more complicated, and some cheap and effective preconditioners cannot be used as the matrix is not available. This cost is even worse in simulations where the system matrix changes every timestep, nonlinear step, or optimization step (in PDE-constrained optimization). Even if we can compute a good preconditioner, we often still see (some) convergence issues, for example, due to PDE coefficients that vary drastically in space, due to strong nonlinearity, highly deformed meshes, and/or poor meshes due to geometric complexity. On the other hand, if we solve the same system many times (time-stepping with a fixed time step) or the system changes slowly, we have substantial opportunity to learn properties of the system from these solves, and we can dynamically update the preconditioner or learn a better preconditioner.

I will discuss strategies to update the preconditioner dynamically using information that is derived during the solve and provide methods for analysis. This is joint work with Kirk Soodhalter and Kapil Ahuja.

Parallel-in-Time for Optimization of PDEs

We will discuss parallel-in-time preconditioning approaches for two different application areas arising from the optimization of PDEs: optimal fluid flow control problems, and mean-field games. The key in each case is to design an iterative algorithm around the best possible choice(s) of fast discrete transforms.

Coffee Break

wherein there will be coffee, tea, and small bites to eat

Stage-Parallel Implicit Runge-Kutta Methods Via Low-Rank Matrix Equation Corrections

Implicit Runge–Kutta (IRK) methods are highly effective for solving stiff ordinary differential equations (ODEs) but can be computationally expensive for largescale problems due to the need of solving coupled algebraic equations at each step. The key idea here is to reformulate a perturbed stage system in a stable way and to retrieve the exact solution through the solution of a Sylvester matrix equation with a known low-rank structure on the right-hand side. We focus on two major IRK families—symmetric and collocation schemes—and extend the methodology to nonlinear settings via a simplified Newton iteration. A set of numerical experiments, including ODEs derived from spatial discretizations of PDEs, confirm the effectiveness of the proposed approach.

Lunch

Lunch catered from Sweet Cecily

Weighted GMRES with preconditioning and deflation

In this talk I will present some new results related to the analysis and acceleration of GMRES. Together with Daniel Szyld (Temple University), we have analyzed the convergence of weighted, preconditioned and deflated GMRES. The emphasis is on positive definite linear systems Ax = b and Hermitian positive definite preconditioners. We show that if a good preconditioner is known for the Hermitian part of A, then some of the properties carry over to the whole system. We also propose a new deflation space. Together with Pierre Matalon, we have shown that any two residual curves can be simultaneously satisfied by left and right preconditioned GMRES. The result is independent of the eigenvalues of the preconditioned operator.

Parallelization of all-at-once preconditioned solvers for time-dependent PDEs

Modern high performance computers provide tremendous compute power by utilizing large amounts of cores. As a consequence, traditional spatial parallelization schemes lead to a saturation of the speedup more often. This motivated the use of parallelization in the time diretion, as well. From a numerical linear algebra viewpoint in this case it is natural to consider all-at-once systems. In many cases a block epsilon-circulant preconditioner can be applied to speed up the convergence of either GMRES or MINRES. In this case multiplication and inversion can be carried out in almost optimal, i.e., O(nlog n), complexity by using the FFT.

While the study of the preconditioners is extensive, the actual parallel implementation is studied very little. One option to implement these kinds of methods is the usage of a parallel FFT. Yet, an efficient parallelization of the FFT is relatively difficult given that it requires a lot of communication in comparison to very few arithmetic operations. This is one reason why multi-dimensional FFTs usually transpose the data such that the individual 1D-FFTs can be carried out sequentially. For the preconditioners considered, in general only 1D-FFTs are needed.

Given the results obtained using efficient multi-dimensional FFTs on parallel computers as an alternative to a direct parallelization of the FFT we propose to transpose the data such that sequential 1D-FFTs can be used and to transpose it back before solving the individual blocks.

We will provide an overview over the different preconditioners that can be implemented in the proposed way, present the parallelization approach in detail, discuss the solution of the blocks and demonstrate the achieved performance.

Closing

Closing remarks

Collaboration Time

Augmented Krylov subspace methods

In many High-Performance Computing applications, it is often of interest to accelerate an iterative Krylov subspace solver, by augmenting the Krylov subspace with some additional subspace. We discuss a new framework [Soodhalter, de Sturler, Kilmer - GAMM Mitteilungen 2020] which describes all augmented Krylov subspace methods in terms of applying a standard Krylov subspace iteration to a projected problem, and taking an additional correction to the solution from the augmentation space. We then discuss new open questions arising from this framework.

Comparison of Haar Wavelet Sparsity and Anisotropic Total Variation in Image Reconstruction for Cone Beam Computed Tomography

Cone beam computed tomography (CBCT) is an increasingly popular imaging modality in which X-ray data of an object is collected using a cone-shaped beam, and a reconstruction algorithm is applied to create a 3D volumetric image. CBCT has a broad range of applications, especially in medical imaging of the head and neck region and of the extremities.

Analytical reconstruction algorithms, while fast and well understood, suffer from noise and cone beam artifacts, and perform poorly with undersampled data, and there is a growing demand for new, more flexible techniques. In this work, we have developed iterative reconstruction algorithms using two different regularization techniques: the Haar wavelet transform and anisotropic total variation. The reconstruction is computed using the primal-dual fixed point (PDFP) algorithm, an effective and flexible optimization technique which allows including additional useful properties such as nonnegative X-ray attenuation. Multiple test cases for the algorithms are provided using both simulated and measured X-ray data.

Learning a microlocal prior for limited-angle tomography

Digital breast tomosynthesis is an application of limited-angle tomography, which is a highly ill-posed inverse problem. Due to the limited-angle imaging geometry, reconstructions suffer from severe stretching of features along the central direction of projections, leading to poor slice separation. We propose a method for learning a boundary estimate for features. This estimate can be presented on top of the reconstruction, indicating the true form and extent of features. Learning the boundary estimate is based on directional edge detection that is implemented using complex wavelets and morphological operations. By using deep learning, we first extract the visible part of the wavefront set and then extend it to the full domain, filling in the parts of the wavefront set that would otherwise be hidden due to the lack of measurement directions. The resulting singular support gives the boundary estimate curve.

Material-separating regularizer formulti-energy X-ray tomography

Dual-energy X-ray tomography is considered in a context where the target under imaging consists of two distinct materials. The materials are assumed to be possibly intertwined in space, but at any given location there is only one material present. Further, two X-ray energies are chosen so that there is a clear difference in the spectral dependence of the attenuation coefficients of the two materials. A novel regularizer is presented for the inverse problem of reconstructing separate tomographic images for the two materials. A combination of two things, (a) non-negativity constraint, and (b) penalty term containing the inner product between the two material images, promotes the presence of at most one material in a given pixel. A preconditioned interior point method is derived for the minimization of the regularization functional. Numerical tests with digital phantoms suggest that the new algorithm outperforms the baseline method, Joint Total Variation regularization, in terms of correctly material-characterized pixels. While the method is tested only in a two-dimensional setting with two materials and two energies, the approach readily generalizes to three dimensions and more materials. The number of materials just needs to match the number of energies used in imaging.

Using persistent homology to improve tomographic reconstructions

The limited-angle tomography has great potential in both medical and industrial imaging. However, the reconstructions have stretching artifacts, and the singularities, i.e., the boundaries of the target object, cannot be detected stably. The goal of our research is to develop a method for detecting unknown singularities from known ones using persistent homology. Thus, this study brings together inverse problems with algebraic topology.

The complex wavelets provide a method for finding stable singularities and dividing them into six subsets based on their directions. The true singularity support forms cycles, but we know only parts of the cycles when using limited angle tomography. Using prior information about directions, it is possible to estimate unknown singularities based on known singularities. Persistent homology can identify cycles when known singularities and estimated singularities occur together. This will reduce our singular support estimation and result in more accurate reconstructions.

Register to attend at THIS FORM so that we know you are coming.

Ireland is not a member of the Schengen Area and a Schengen visa is not valid for travel to Ireland. A passport is required for travel to Ireland from any other state (including the United Kingdom where it is not an immigration requirement, but airlines require identification documents and many will only accept a passport).

Please see this advice from the Irish government on whether or not you need to apply for a visa in advance of travel to Ireland to attend an academic conference, according to your nationality as shown on your passport or travel document. Participants who do need a visa can find advice there on the online application process and the likely duration of the wait for approval.

Dublin Airport is the main arrival point for international visitors to Ireland and has daily direct flights to many locations worldwide, particularly in Europe and North America. It is connected by multiple daily flights to the major European hubs such as Frankfurt, Heathrow, Amsterdam and Charles de Gaulle.

Dublin Airport provides assistance for travellers with constrained mobility, which can be arranged in advance.

Travellers from Dublin Airport to the United States clear US Immigration and Customs in Dublin, which requires a little extra time in Dublin Airport before scheduled departure.

There are frequent buses which depart from Dublin airport to the center of the city. More information about the different bus services available from the airport can be found at the airport website. One can also catch a taxi from the taxi ranks. Please follow signs in the airport to an official taxi rank booth.

The entire workshop will take place in the Hamilton Building, which is located on the far east side of Trinity's main campus. Although the Hamilton Buiding has an address on Westland Row, it cannot be entered by the public from this side. Rather, one should enter campus either at the Pearse Street entrance or the Lincoln Gate entrance. From either entrance, one can enter the Hamilton Building at one of two points (here or here).

All talks will take place in the Salmon lecture theater. Salmon is on the Ground Floor roughly halfway between the two ground-floor entrances mentioned in the previous paragraph.