Registration Room: Sala degli Stemmi, Scuola Normale Superiore

Opening Remarks Room: Sala degli Stemmi, Scuola Normale Superiore

Why Parallel in Time (PinT) methods are different for Parabolic and Hyperbolic Problems

• Martin J. Gander (University of Geneva)

Room: Sala degli Stemmi, Scuola Normale Superiore Space-time parallel methods, also known more recently under the name PinT methods, have received a lot of attention over the past two decades, driven by the parallel hardware architectures that have now millions of cores, and acceleration often saturates when one parallelizes in space only. Research focus has therefore shifted to trying to parallelize also the time direction. However, for parallelization, time is very different from space, because evolution problems satisfy a causality principle: the future is dependent on the past, and not the other way round. It is therefore not clear a priori if useful numerical work can be done simultaneously in the near and far future. I will first show in my presentation why intuitively there is a fundamental difference when one parallelizes hyperbolic or parabolic problems in space-time, and which key properties need to be taken into account to be successful. I will then give examples of PinT methods which use the physical properties of the evolution problem to their advantage to parallelize in space-time. For hyperbolic problems, effective PinT methods are Domain Decomposition methods of Waveform Relaxation type, culminating in Unmapped Tent Pitching methods. Very powerfull techniques are also the ParaDiag methods, and direct time parallel methods like ParaExp. Several of these methods can also be very effectively used for parabolic problems, but for such problems there are also very successful multilevel methods, like Parareal, and the currently best ones are space-time multigrid methods. These multilevel methods struggle however when applied to hyperbolic problems.

Parallel-In-Time solver for the all-at-once Runge--Kutta discretization

• Luca Bergamaschi (University of Padua)

Room: Sala degli Stemmi, Scuola Normale Superiore Time-dependent PDEs arise quite often in many scientific areas, such as mechanics, biology, economics, or chemistry. Of late, researchers have devoted their effort in devising parallel-in-time methods for the numerical solution of such PDEs, adding a new dimension of parallelism and allowing to speed-up the solution process on modern supercomputers. In this talk, we present a fully parallelizable preconditioner for the all-at-once linear system arising when employing a Runge-Kutta (RK) method in time. The resulting system is solved iteratively for the numerical solution and for the stages. The proposed preconditioner results in a block-diagonal solve for the stages, accelerated by a novel block-preconditioner based on the SVD of the RK coefficient matrix, and a Schur complement obtained by solving again systems for the stages. Parallel results show the robustness of the preconditioner with respect to the discretization parameters and to the number of stages, as well as very promising scalability and parallel efficiency indices.

A new ParaDiag time-parallel time integration method

• Davide Palitta (University of Bologna)

Room: Sala degli Stemmi, Scuola Normale Superiore Many parallel-in-time (PinT) paradigms have been developed in the last decades to efficiently solve time-dependent partial differential equations (PDEs). In this talk we focus on the ParaDiag scheme whose most peculiar feature consists in the explicit diagonalization of the matrix stemming from the adopted time integrator. However, certain classes of time integrators lead to a discrete operator which is not diagonalizable. For instance, this is the case of Backward Differentiation Formulas (BDFs) like the backward Euler scheme. With the aim of overcoming such a drawback, different approaches have been developed in the literature. In this talk, we illustrate a novel technique. In particular, we show how to exploit the possible circulant-plus-low-rank structure of the discrete time integrator to design a new, successful implementation of the ParaDiag paradigm. Notice that this peculiar structure arises in many different families of time integrators. The efficiency of our original scheme is displayed by several numerical examples.

A preconditioned MINRES method for block lower triangular Toeplitz systems with application to solving all-at-once systems of linear time-dependent PDEs

• Sean Hon (Hong Kong Baptist University)

Room: Sala degli Stemmi, Scuola Normale Superiore In this study, a novel preconditioner based on the absolute-value block α-circulant matrix approximation is developed, specifically designed for nonsymmetric dense block lower triangular Toeplitz (BLTT) systems that emerge from the numerical discretization of evolutionary equations. Our preconditioner is constructed by taking the absolute value of a block α-circulant matrix approximation to the BLTT matrix. To apply our preconditioner, the original BLTT linear system is converted into a symmetric form by applying a time-reversing permutation transformation. Then, with our preconditioner, the preconditioned minimal residual method (MINRES) solver is employed to solve the symmetrized linear system. With a properly chosen α, the eigenvalues of the preconditioned matrix are proven to be clustered around ±1 without any significant outliers. With the clustered spectrum, we demonstrate that the preconditioned MINRES solver for the preconditioned system has a convergence rate independent of system size. Our preconditioner can be implemented in a parallel-in-time manner. The efficacy of the proposed preconditioner is corroborated by our numerical experiments, which reveal that it attains optimal convergence.

Global space-time low-rank methods for the time-dependent Schrödinger equations

• Virginie Ehrlacher (Ecole des Ponts ParisTech)

Room: Sala degli Stemmi, Scuola Normale Superiore The aim of this talk is to present novel global space-time methods for the approximation of the time-dependent Schrödinger equation, using Kato theory. The latter can be used in conjunction with low-rank tensor formats (such as Tensor Trains for instance) to derive new variational principles to compute dynamical low-rank approximations of the solution, which are different from the Dirac-Frenkel principle. One significant advantage of this new variational formulation is that the existence of a dynamical low-rank approximation for any finite-time horizon can be proved, whereas dynamical low-rank approximations constructed with the Dirac-Frenkel principle can usually be porved to exist only locally in time. Illustrative numerical results will be presented to highlight the differences between the dynamical low-rank approximations obtained with these different approaches. This is joint work with Clément Guillot and Mi-Song Dupuy.

Parallel time-dependent variational principle algorithm for matrix product states

• Sergey Dolgov (University of Bath)

Room: Sala degli Stemmi, Scuola Normale Superiore Combining the time-dependent variational principle (TDVP) algorithm with the parallelization scheme introduced by Stoudenmire and White for the density matrix renormalization group (DMRG), we present the first parallel matrix product state (MPS) algorithm capable of time evolving one-dimensional (1D) quantum lattice systems with long-range interactions. We benchmark the accuracy and performance of the algorithm by simulating quenches in the long-range Ising and XY models. We show that our code scales well up to 32 processes, with parallel efficiencies as high as 86%.

Parallel implementation of block circulant type preconditioner for all-at-once systems of linear time-dependent PDEs

• Ryo Yoda (The University of Tokyo)

Room: Sala degli Stemmi, Scuola Normale Superiore Parallel-in-time approaches solve all-at-once systems obtained by solving all time-dependent PDEs at once in order to extract temporal parallelism. Assuming linear and constant-in-time integrators for all time steps, the resulting systems have a block Toeplitz structure. Consequently block circulant preconditioners have attracted much attention for these systems. In particular, block epsilon-circulant preconditioners, introducing a weighting coefficient, have achieved convergence independent from the spatial size, thus they are promising parallel-in-time approaches. This work focuses on parallel implementations of block circulant type preconditioners. The primary operations of these preconditioners are FFTs on one-dimensional time-step-sized vectors and solving spatial-sized linear systems with complex-valued coefficient matrices. We use FFTW and Trilinos packages with their MPI implementation and investigate their parallel performance. Additionally, we propose an alternative parallelization strategy to execute the one-dimensional FFTs. Numerical experiments demonstrate good scaling behavior for linear diffusion and advection-diffusion problems.

Efficiently solving non-linear time-dependent problems from the geosciences

• Gabriel Wittum (King Abdullah University of Science and Technology)

Room: Sala degli Stemmi, Scuola Normale Superiore The presentation investigates the efficient use of a linearly implicit stiff integrator for the numerical solution of density driven flow problems. Upon choosing a one-step method of extrapolation type (code LIMEX), the use of full Jacobians and reduced approximations are discussed. Numerical experiments include nonlinear density flow problems such as diffusion from a salt dome (2D), a (modified) Elder problem (3D), the saltpool benchmark (3D) and a real life salt dome problem (2D). The arising linear equations are solved using either a multigrid preconditioner from the software package UG4 or the sparse matrix solver SuperLU. Based on these component, this work devises guidelines for the design of an efficient solver.

Regularized dynamical parametric approximation

• Christian Lubich (University of Tuebingen)

Room: Sala degli Stemmi, Scuola Normale Superiore This paper studies the numerical approximation of solutions to initial value problems of high-dimensional ordinary differential equations or evolutionary partial differential equations such as the Schr\"odinger equation by nonlinear parametrizations $u(t)=\Phi(q(t))$ with time-dependent parameters $q(t)$, which are to be determined in the computation. The motivation comes from approximations in quantum dynamics by multiple Gaussians and approximations of various dynamical problems by tensor networks and by neural networks. In all these cases, the parametrization is typically irregular: the derivative $\Phi'(q)$ can have arbitrarily small singular values and may have varying rank. We derive approximation results for a regularized approach in the time-continuous case as well as in time-discretized cases. With a suitable choice of the regularization parameter and the time stepsize, this approach can still be successfully applied in such irregular situations, even if it runs counter to the basic principle in numerical analysis to avoid solving ill-posed subproblems when aiming for a stable algorithm. Numerical experiments with sums of Gaussians for approximating laser-induced quantum dynamics and with neural networks for approximating the flow map of a system of ordinary differential equations illustrate and complement the theoretical results. The talk is based on joint work with Caroline Lasser, Jörg Nick and Michael Feischl.

Randomized low-rank approximation for time- and parameter-dependent problems

• Daniel Kressner (EPFL)

Room: Sala degli Stemmi, Scuola Normale Superiore A variety of applications gives rise to matrices $A(t)$ that admit good low-rank approximation for each time (or parameter value) $t$, including dynamical systems, spectral density estimation, and Gaussian process regression, In this talk, we discuss the benefits of randomized methods for approximating $A(t)$ simultaneously for many values of $t$. We describe and analyze parameter-dependent extensions of two popular randomized algorithms, the randomized singular value decomposition and the generalized Nyström method. Both, the theoretical results and numerical experiments, show that these methods reliably return quasi-best low-rank approximations. This talk is based on joint work with Hysan Lam.

Adaptive rational Krylov methods for exponential Runge--Kutta integrators

• Kai Bergermann (TU Chemnitz)

Room: Sala degli Stemmi, Scuola Normale Superiore We consider the solution of large stiff systems of ordinary differential equations with explicit exponential Runge--Kutta integrators. These problems arise from semi-discretized semi-linear parabolic partial differential equations on continuous domains or on inherently discrete graph domains. A series of results reduces the requirement of computing linear combinations of $\varphi$-functions in exponential integrators to the approximation of the action of a smaller number of matrix exponentials on certain vectors. State-of-the-art computational methods use polynomial Krylov subspaces of adaptive size for this task. They have the drawback that the required number of Krylov subspace iterations to obtain a desired tolerance increase drastically with the spectral radius of the discrete linear differential operator, e.g., the problem size. We present an approach that leverages rational Krylov subspace methods promising superior approximation qualities. We prove a novel a-posteriori error estimate of rational Krylov approximations to the action of the matrix exponential on vectors for single time points, which allows for an adaptive approach similar to existing polynomial Krylov techniques. We discuss pole selection and the efficient solution of the arising sequences of shifted linear systems by direct and preconditioned iterative solvers. Numerical experiments show that our method outperforms the state of the art for sufficiently large spectral radii of the discrete linear differential operators. The key to this are approximately constant numbers of rational Krylov iterations, which enable a near-linear scaling of the runtime with respect to the problem size. We will also discuss the special algebraic structure of the supra-Laplacian operator appearing in multilayer network analysis.

Energy-preserving splitting integration for Hamiltonian Monte Carlo method with adaptive tuning

• Fasma Diele (IAC)

Room: Sala degli Stemmi, Scuola Normale Superiore Splitting schemes provide a promising alternative to the classical Stormer-Verlet method in Hamiltonian Monte Carlo (HMC) methodology. Within the family of one-parameter second-order splitting procedures, we demonstrate that using a designated function of the free parameter to select the step size ensures stability and Hamiltonian preservation when sampling from Gaussian distributions. This guarantees no sample rejections in the HMC process, a key factor for superior performance compared to recent similar methods. The effectiveness of the proposed approach for sampling from general non-Gaussian distributions is assessed, incorporating a simple adaptive selection technique for the free parameter to improve HMC performance. Benchmark examples from literature and experiments, including the Log-Gaussian Cox process and Bayesian Logistic Regression, highlight the effectiveness of the approach.

Detecting the numerical of ill posedness in delay differential equations

• Miryam Gnazzo (GSSI)

Room: Sala degli Stemmi, Scuola Normale Superiore Given a set of matrices $A_i \in \mathbb{C}^{n \times n}$ and a set of analytic functions $f_i : \mathbb{C} \mapsto \mathbb{C}$, we consider a regular matrix-valued function $\mathcal{F}(\lambda)=\sum_{i=0}^d f_i(\lambda) A_i$, that is $\det \left( \mathcal{F}(\lambda) \right)$ is not identically zero for $\lambda \in \mathbb{C}$. An interesting problem consists in the computation of the nearest singular function $\widetilde{\mathcal{F}}(\lambda)=\sum_{i=0}^d f_i(\lambda) \left( A_i + \Delta A_i \right)$, with respect to the Frobenius norm. For example, this problem has particular importance in the context of delay differential algebraic equations, where a function in the form $\mathcal{D}(\lambda)=\lambda E - A - B e^{-\tau \lambda}$ is studied. Indeed in this setting and in presence of small delays $\tau$, the ill posedness of the problem may be connected with the numerical singularity of the function $\mathcal{D}(\lambda)$, even if the pencil $\lambda E -A$ is regular. We will provide a general overview of the problem, describing the possible issues connected with the lack of robustness of the differential equation, associated with destabilizing perturbations of $\mathcal{D}(\lambda)$. Moreover we propose a method for the numerical approximation of the function $\mathcal{F}(\lambda)$, which rephrases the matrix nearness problem for the matrix-valued function into an equivalent optimization problem. Nevertheless this problem turns out to be highly non-convex. To solve it, we propose a two level procedure, which introduces a constrained gradient system of differential equations in the inner iteration and a Newton-like method for the optimization of the perturbation size in the outer one. This is a joint work with Nicola Guglielmi (GSSI).

High-order conservative and accurately dissipative numerical integrators via finite elements in time

• Patrick Farrell (University of Oxford)

Room: Sala degli Stemmi, Scuola Normale Superiore Numerical methods for the simulation of transient systems with structure-preserving properties are known to exhibit greater accuracy and physical reliability, in particular over long durations. These schemes are often built on powerful geometric ideas for broad classes of problems, such as Hamiltonian or reversible systems. However, there remain difficulties in devising higher-order-in-time structure-preserving discretizations for nonlinear problems, and in conserving higher-order invariants. In this work we propose a general framework for the construction of structure-preserving timesteppers via finite elements in time and the systematic introduction of auxiliary variables. The framework reduces to Gauss methods where those are structure-preserving, but extends to generate arbitrary-order structure-preserving schemes for nonlinear problems, and allows for the construction of schemes that conserve multiple higher-order invariants. We demonstrate the ideas by devising novel schemes that exactly conserve all known invariants of the Kepler and Kovalevskaya problems, high-order energy-conserving schemes for the compressible Navier--Stokes equations, and high-order energy- and helicity-preserving schemes for the incompressible magnetohydrodynamics equations.

Convolution Quadrature for the quasilinear subdiffusion equation

• Maria Lopez-Fernandez (University of Malaga)

Room: Sala degli Stemmi, Scuola Normale Superiore We construct a Convolution Quadrature (CQ) scheme for the quasilinear subdiffusion equation and supply it with the fast and oblivious implementation. In particular we find a condition for the CQ to be admissible and discretize the spatial part of the equation with the Finite Element Method. We prove the unconditional stability and convergence of the scheme and find a bound on the error. As a passing result, we also obtain a discrete Grönwall inequality for the CQ, which is a crucial ingredient of our convergence proof based on the energy method. The paper is concluded with numerical examples verifying convergence and computation time reduction when using fast and oblivious quadrature.

High-order conservative and accurately dissipative numerical integrators via finite elements in time

• Boris Andrews (University of Oxford)

Room: Sala degli Stemmi, Scuola Normale Superiore Numerical methods for the simulation of integrable systems with conservative properties are known to exhibit greater accuracy and physical reliability, in particular over long durations. One important class of schemes is that of Gauss methods, a class of symplectic Runge–Kutta methods. In this talk we propose an alternative, general framework for the construction of conservative schemes via finite elements in time and the systematic introduction of auxiliary variables. For linear problems with quadratic invariants, the scheme is provably equivalent to Gauss methods. However, the alternative framework extends to nonlinear systems with potentially multiple non-quadratic invariants. The framework also allows for the construction of numerical methods that accurately preserve dissipation structures, e.g. energy dissipation in the incompressible Navier–Stokes equations. We demonstrate the ideas by devising novel schemes that exactly conserve all known invariants of Hamiltonian systems with high-order invariants, and mass– and energy–conserving schemes for the compressible Navier–Stokes equations. This talk complements and extends the plenary talk given by Patrick Farrell.

Geometric integration meets data-driven dynamical systems

• Matthew Colbrook (University of Cambridge)

Room: Sala degli Stemmi, Scuola Normale Superiore Koopman operators globally linearise nonlinear dynamical systems, and their spectral information can be a powerful tool for analysing and decomposing such systems. There has been a recent flurry of activity in data-driven computations of the spectral properties of Koopman operators. However, Koopman operators are infinite-dimensional, making the computation of their spectral information a considerable challenge. In this talk, we will combine ideas from geometric integration with dynamic mode decomposition (DMD), one of the most popular algorithms for analysing Koopman operators. We introduce \textit{measure-preserving extended dynamic mode decomposition} (\texttt{mpEDMD}), the first Galerkin method whose eigendecomposition converges to the spectral quantities of Koopman operators for general measure-preserving dynamical systems. \texttt{mpEDMD} is a data-driven algorithm based on an orthogonal Procrustes problem that enforces measure-preserving truncations of Koopman operators using a general dictionary of observables. It is flexible and easy to use with any pre-existing DMD-type method and with different data types. Enforcing this structure is crucial to its convergence and qualitative behaviour. We prove the convergence of \texttt{mpEDMD} for projection-valued and scalar-valued spectral measures, spectra, and Koopman mode decompositions. For the case of delay embedding (Krylov subspaces), our results include the first convergence rates of the approximation of spectral measures as the size of the dictionary increases. We demonstrate \texttt{mpEDMD} on a range of challenging examples, its increased robustness to noise compared with other DMD-type methods, and its ability to capture the energy conservation and cascade of a turbulent boundary layer flow with a Reynolds number $> 6 \times 10^4$ and a state-space dimension $> 10^5$.

Poster Session and Apero Room: Dipartimento di Matematica, main hall Pietro Deidda Arturo De Marinis Alejandro Escorihuela-Tomàs Sebastian Götschel Stefano Sicilia

A dynamical systems view to deep learning: contractivity and structure preservation

• Elena Celledoni (NTNU)

Room: Sala degli Stemmi, Scuola Normale Superiore The (discrete) optimal control point of view to neural networks offers an interpretation of deep learning from a dynamical systems and numerical analysis perspective and opens the way to mathematical insight. In this talk we discuss topics of structure preservation and their use in deep learning. Some deep neural networks can be designed to have desirable properties such as invertibility and group equivariance or can be adapted to problems of manifold value data. Contractivity is identified as a desirable property for stability and robustness of neural networks. We discuss classical results of contractivity of numerical ODE integrators, applications to neural networks and recent extensions to Riemannian manifolds.

Variational (dynamical) discretisations of metriplectic systems

• Damiano Lombardi (INRIA)

Room: Sala degli Stemmi, Scuola Normale Superiore In this work we present a variational structure preserving discretisation method for metriplectic systems. These are composed of two parts: a Hamiltonian flow and a dissipation. The dissipation is built based on the Casimirs of the Hamiltonian system. A weak formulation of the Hamiltonian flow is defined. Based on this, a variational semi-discretisation in space is performed, leading to a finite dimensional Hamiltonian system. We discuss the possibility of introducing a dynamical variational discretisation by leveraging the geometrical description of the weak Hamiltonian flow, in the spirit of the Dynamical Low Rank method. Work in collaboration with Cecilia Pagliantini.

A new fast numerical method for the generalized Rosen-Zener model

• Stefano Pozza (Charles University)

Room: Sala degli Stemmi, Scuola Normale Superiore In quantum mechanics, the Rosen-Zener model represents a two-level quantum system. Its generalization to multiple degenerate sets of states leads to larger non-autonomous linear systems of ordinary differential equations (ODEs). We propose a new method for computing the solution operator of this system of ODEs. This new method is based on a recently introduced expression of the solution in terms of an infinite matrix equation, which can be efficiently approximated by combining truncation, fixed point iterations, and low-rank approximation. This expression is possible thanks to the so-called ∗-product approach for linear ODEs. In the numerical experiments, the new method’s computing time scales linearly with the model’s size. We provide a first partial explanation of this linear behavior. Joint work with Christian Bonhomme and Niel Van Buggenhout.

Lie-Poisson discretization for incompressible magnetohydrodynamics on the sphere

• Michael Roop (Chalmers University of Technology)

Room: Sala degli Stemmi, Scuola Normale Superiore We give a structure preserving spatio-temporal discretization for incompressible magnetohydrodynamics (MHD) on the sphere. Discretization in space is based on the theory of geometric quantization, which yields a spatially discretized analogue of the MHD equations as a finite-dimensional Lie-Poisson system on the dual of the magnetic extension Lie algebra, for which we develop structure preserving time discretizations. The full method preserves the underlying geometry, namely the Lie-Poisson structure and all the Casimirs. To showcase the method, we apply it to two models for magnetized fluids: incompressible magnetohydrodynamics and Hazeltine's model. This is a joint work with Klas Modin, arXiv:2311.16045.

Variational discretizations of ideal magnetohydrodynamics in smooth regime using finite element exterior calculus

• Valentin Carlier (MPI)

Room: Sala degli Stemmi, Scuola Normale Superiore We propose a new class of finite element approximations to ideal compressible magnetohydrodynamics equations in smooth regime. Our discretizations are built via a discrete variational principle mimicking the continuous Euler-Poincaré principle, with vector fields represented by their action as Lie derivatives on differential forms, to further exploit the geometrical structure of the problem. The resulting semi-discrete approximations are shown to conserve the total mass, entropy and energy of the solutions. In addition the divergence-free nature of the magnetic field is preserved in a pointwise sense, and the scheme is reversible at the fully discrete level. Numerical simulations are conducted to verify the accuracy of our approach and its ability to preserve the semi-discrete invariants for several test problems. An eponym paper will be uploaded to arxiv in the coming days.

Concluding remarks Room: Sala degli Stemmi, Scuola Normale Superiore