ABSTRACT. In the literature, the convergence analysis of multigrid methods for structured matrix-sequences has been obtained in a compact and elegant form in the case of sequences associated with a (possibly multivariate) scalar-valued function $f$ [1,6]. However, in the case where the entries are generic $s\times s$ matrices instead of scalars the general interest mainly regards specific applications, that is, the focus concerns algorithmic proposals with very marginal attention to general theoretical results. Hence, there is still a substantial lack of an effective strategy for choosing the grid transfer operators and of a rigorous convergence analysis.

The first aim of the talk is to generalize the convergence results of multigrid algorithms applied on linear systems where the coefficient matrix sequence is a block circulant or a block Toeplitz matrix sequence associated with a matrix-valued symbol $\mathbf{f}$ [3]. First, we focus on the crucial choice of conditions that should be fulfilled by the trigonometric polynomial $\mathbf{p}$ used to construct the projector in order to ensure the optimality of the two-grid method [2,3].

In particular, the latter are expressed in terms of eigenvectors associated with the ill-conditioned subspace and permit to extend the class of suitable trigonometric polynomials that can be chosen to construct the projectors.

Then, we generalize the convergence analysis to the V-cycle method proving a linear convergence rate under more restrictive conditions on $\mathbf{p}$, which resemble those given in the scalar case.

Moreover, we prove the efficiency of some standard projectors [4,5] whose associated polynomial, indeed, fulfills our derived conditions. Finally, using a tensor product argument, we provide a strategy to construct efficient V-cycle procedures in the most general block multilevel setting.

[1] A. Aricò, M. Donatelli. A V-cycle multigrid for multilevel matrix algebras: proof of optimality. Numer. Math., 105 (4) (2007), 511–547. [2] M. Bolten, M. Donatelli, P. Ferrari, I. Furci. A symbol based analysis for multigrid methods for Block-Circulant and Block-Toeplitz Systems. Under review. [3] M. Donatelli, P. Ferrari, I. Furci, S. Serra–Capizzano, D. Sesana. Multigrid methods for Block-Toeplitz linear systems: Convergence Analysis and Applications. Numer. Linear Algebra Appl., (2021), e2356. [4] P. Ferrari, R. I. Rahla, C. Tablino Possio, S. Belhaj, S. Serra–Capizzano. Multigrid for Q k Finite Element Matrices using a (block) Toeplitz symbol approach. Mathematics, 8 (1-5) (2020). [5] W. Hackbusch. Multigrid methods and applications. vol. 4 of Springer Series in Computational Mathematics Springer-Verlag, Berlin (1985). [6] J. W. Ruge, K. Stüben. Algebraic multigrid. Multigrid Methods, S. McCormick, ed., Frontiers Appl. Math. 3, SIAM, Philadelphia (1987), 73–130.

ABSTRACT. In [1], the authors proposed a general two-grid convergence analysis, proving an optimal convergence rate independent of the matrix size, in the case of positive definite block-Toeplitz matrices with generic blocks. The aim of this talk is to show how to exploit the theoretical findings in [1] to develop and analyse multigrid strategies for the solution of linear systems stemming from rectangular Lagrangian Finite Elements approximations of elliptic partial differential equations.

Firstly, we propose a classical multigrid strategy that follows a functional approach, that is, we define the prolongation operator as the inclusion operator between the coarser and finer functional spaces [2]. We analyse the prolongation matrix as a cut block-Toeplitz matrix and we prove that its symbol satisfies the hypotheses for the two-grid and V-cycle convergence.

Moreover, we present a class of grid transfer operators defined according to the analysis in [1]. Indeed, we explain how to choose the trigonometric polynomial that generates the block-Toeplitz matrix used in the construction of the grid transfer operator focusing only on algebraic considerations on the symbol of the linear system matrix-sequence. Even though we focus on the rectangular Lagrangian Finite Elements stiffness matrices, the presented procedure has a wider interest, since it might be applied to every matrix-sequence that falls into the theoretical setting.

Finally, results of numerical experiments that test all the considered methods are presented, both in one dimension and in higher dimension, showing an optimal behaviour in terms of the dependency on the matrix size and a robustness with respect to the dimensionality.

[1] M. Donatelli, P. Ferrari, I. Furci, S. Serra–Capizzano, D. Sesana. Multigrid methods for Block-Toeplitz linear systems: Convergence Analysis and Applications. Numer. Linear Algebra Appl., (2021), e2356. [2] P. Ferrari, R. I. Rahla, C. Tablino Possio, S. Belhaj, S. Serra–Capizzano. Multigrid for Q k Finite Element Matrices using a (block) Toeplitz symbol approach. Mathematics, 8 (1-5) (2020).

ABSTRACT. In numerical simulations, memory bandwidth often becomes a major performance bottleneck. This is particularly true for linear systems in extreme scale multigrid applications when storing the system matrix may not even be feasible. This has led to the development of matrix-free methods, where the matrix coefficients are computed on the fly. Of course, such methods usually result in computational overhead for executing the matrix assembly repeatedly.

It can be shown that, for suitable finite element meshes, the stiffness matrix can be approximated by polynomial surrogate functions. The efficiency of this method relies on balancing the accuracy of the surrogate polynomials and the computational cost incurred, including memory accesses. This work presents a concise analysis of the potential speedup that can be achieved. Numerical studies with the Poisson- and Stokes equations support the theoretical results.

Both, the analysis and the numerical evidence confirm the high potential of the method.

ABSTRACT. We focus on a two-dimensional time-space diffusion equation with fractional derivatives in space. The use of Crank-Nicolson in time and finite differences in space leads to dense Toeplitz-like linear systems. Multigrid strategies that exploit such structure are particularly effective when the fractional orders are both close to 2. We seek to investigate how structure-based multigrid approaches can be efficiently extended to the case where only one of the two fractional orders is close to 2, i.e., when the fractional equation shows an intrinsic anisotropy. Precisely, we design a multigrid (block-banded–banded-block) preconditioner whose grid transfer operator is obtained with a semi-coarsening technique and that has relaxed Jacobi as smoother. The Jacobi relaxation parameter is estimated by using an automatic symbol-based procedure. A further improvement in the robustness of the proposed multigrid method is attained using the V-cycle with semi-coarsening as smoother inside an outer full-coarsening. Several numerical results confirm that the resulting multigrid preconditioner is computationally effective and outperforms current state of the art techniques.

ABSTRACT. We present a p-multigrid accelerated pseudospectral scheme for efficient and scalable simulation of dispersive and steep nonlinear water waves in a numerical wave tank setup based on potential flow theory. The numerical scheme is based on a Fourier collocation method for the discretization of the free surface equations, and a hybrid Fourier-Chebyshev collocation method for the discretization of the Laplace equation in the fluid. To help ensure stability of the scheme, the free surface equations have been de-aliased by employing Orszag’s classical rule. The computational bottleneck of the problem is the Laplace equation which is solved iteratively using a p-multigrid accelerated defect correct method. The acceleration is done through pre-conditioning using a p-multigrid method for the linearised wave-problem. Through numerical analysis we examine the efficiency and scalability of the proposed scheme and give examples of numerical bench- marks for nonlinear water waves. The work is linked to ongoing efforts towards establishing improved robust tools for large-scale wave-simulations for offshore engineering applications.

ABSTRACT. Unfitted finite element methods are quite popular for modeling the complex domains, as they allow us to use relatively simple meshes, e.g., structured meshes, for the discretization of the geometries. In these methods, the finite element (FE) spaces associated with the background meshes are modified to capture the information about the geometry. We consider the discretization of contact problems in the unfitted FE framework and employ the Lagrange multiplier method to enforce the contact condition on unfitted boundaries or interfaces. In order to achieve a stable Lagrange multiplier space, we follow the vital-vertex strategy and create a multiplier space that satisfies the discrete inf-sup condition. The benefit of the Lagrange multiplier formulation over, e.g., Nitsche's method, is that the discretization remains unchanged for either equality or inequality constraints and only the solution method has to be modified. However, we need to solve a linearly constrained convex minimization problem instead of a convex unconstrained minimization problem.

In order to solve this linearly constrained minimization problem, we propose a tailored generalized multigrid method. The proposed multigrid method is an extension of the monotone multigrid method, for solving more generic linear bound constraints. For tackling the constraints, we introduce a technique to decouple the linear constraints by projecting them onto a new basis. The decoupled constraints are then handled by a modified projected Gauss-Seidel method, which will be used as a smoother within our MG method. The modified projected Gauss-Seidel methods can enforce the linear constraints locally and monotonically minimize the energy functional, which ensures the global convergence of our method. For the transfer operators, we use a so-called pseudo-$L^2$-projection in the unfitted finite element framework. We shall demonstrate the robustness and the efficiency of the proposed multigrid method in a series of numerical examples and show its level independent convergence property.

ABSTRACT. We develop an augmented Lagrangian preconditioner for a discretisation based on the Scott–Vogelius finite element pair for the velocity and pressure for a system describing the thermal convection of an incompressible fluid with non-Newtonian rheology. The preconditioner involves a specialised multigrid algorithm that makes use of a space decomposition that captures the kernel of the divergence and non-standard intergrid transfer operators. The preconditioner exhibits robust convergence behaviour, even for systems including viscous dissipation and temperature-dependent viscosity and heat conductivity.

ABSTRACT. Low-order finite element linear systems can be used to construct effective preconditioners to the linear systems arising from higher-order discretizations of the Stokes equations. In this paper, we investigate the use of geometric multigrid based on the Q1isoQ2/Q1 discretization of the Stokes operator as a preconditioner to a better-known Q2/Q1 discretization of the Stokes system. We utilize local Fourier analysis to optimize the damping parameters for the Vanka and Braess-Sarazin relaxation schemes and to achieve robust convergence. These results are then verified and compared against the measured multigrid performance. While in principle geometric multigrid can be applied directly to the Q2/Q1 system, our ultimate motivation is to apply algebraic multigrid within solvers for Q2/Q1 systems, which will be considered in a companion paper.

ABSTRACT. Multigrid methods are one of the most efficient techniques for solving linear systems arising from Partial Differential Equations (PDEs) and graph Laplacians from machine learning applications. One of the key components of multigrid is smoothing, which aims at reducing high-frequency errors on each grid level. However, finding optimal smoothing algorithms is problem-dependent and can impose challenges for many problems. In this paper, we propose an efficient adaptive framework for learning optimized smoothers from operator stencils in the form of convolutional neural networks (CNNs). The CNNs are trained on small-scale problems from a given type of PDEs based on a supervised loss function derived from multigrid convergence theories, and can be applied to large-scale problems of the same class of PDEs. Numerical results on anisotropic rotated Laplacian problems demonstrate improved convergence rates and solution time compared with classical hand-crafted relaxation methods.

ABSTRACT. The convergence rate of classical domain decomposition methods in general deteriorates severly for large variations or jumps within material properties, resulting in large variations in the spectrum of the system. To retain the robustness for such highly heterogeneous problems, the coarse space can be enriched by additional coarse basis functions, which are computed from the solutions of local generalized eigenvalue problems. However, the set-up and the solution of the generalized eigenvalue problems typically takes up a significant part of the total time to solution, especially in three dimensions. In particular, for many realistic model problems, only the solution of the eigenvalue problems on a small part of the interface is necessary to design a robust coarse space. In general, it is difficult to predict a priori which of the eigenvalue problems are needed. Here, we use machine learning techniques to predict where eigenvalue problems have to be solved, often reducing its number significantly. In this talk, we extend our results for two dimensions to three dimensions which requires additional and very specific effort in the generation and preparation of the training data. We provide numerical results for linear diffusion and elasticity problems considering both, regular and graph partitioned decompositions.

ABSTRACT. Deep residual networks (ResNets) are widely used neural network architectures. Their popularity originates from their remarkable performance in complex statistical learning tasks, such as image segmentation or classification. The training of ResNets is typically performed using variants of the stochastic gradient (SGD) method. Although the SGD method has a low computational cost per iteration, it often exhibits poor convergence properties, which rely on the choice of hyper-parameters. In particular, it is important to carefully select a sequence of diminishing step-sizes to ensure convergence to a solution. In this work, we propose to train ResNets using a nonlinear multilevel minimization method, namely the recursive multilevel trust region (RMTR) method. The RMTR method combines the global convergence property of the trust region method and the optimality of the multilevel method. The use of the trust-region method also reduces the dependency on the hyper-parameters, as the sequence of the step-sizes is determined automatically by the trust-region radius. The multilevel hierarchy is constructed by exploiting the fact that the ResNet architecture can be reformulated as a forward Euler discretization of a nonlinear initial value problem (IVP). This allows us to construct the multilevel hierarchy of ResNets with different depths by discretizing the original IVP with larger time-steps. The devised RMTR method operates in a hybrid stochastic-deterministic regime and adaptively adjusts sample sizes during the training process. In this work, we will analyze the convergence behavior of our proposed multilevel minimization method and show its robustness.We will demonstrate the significant speedup achieved by the RMTR method in comparison with the standard SGD method or single-level trust-region method using several numerical examples.

ABSTRACT. Physics-informed neural networks (PINNs) have emerged as a promising tool for solving inverse problems for PDEs. In PINNs, the solution to a PDE is approximated by a deep neural network. Unknown parameters in the PDE and parameters in the neural network are found simultaneously by minimizing a residual, typically via gradient-based optimizers. Due to the high dimensionality and lack of convexity of this optimization problem, it is often difficult to find suitable solutions. In this work, we develop a multi-scale training technique that can improve training. Our method creates a hierarchy of residuals at increasing resolution and alternates the optimizer among the residuals in a multi-grid cycle. We demonstrate that this training method results in better training accuracy for PINNs and control volume PINNs, a PINNs variant designed for conservation laws.

ABSTRACT. When designing multigrid methods for hybridized finite element discretizations, natural injection operators do not exist since the refined mesh contains edges which are not refinements of coarse edges. There, a continuation of the coarse grid function on the skeleton must be defined. We discuss conditions on the stability of such injection operators and propose several variants. Employing these, we obtain a uniformly convergent multigrid algorithm for hybridized DG and mixed finite element methods which uses the same discretization on all levels. We discuss some key points of the analysis and show numerical results to prove feasibility of the algorithm.

ABSTRACT. We consider a saddle point least squares method for discretizing the singularly perturbed reaction-diffusion problem. Choices for discrete stable spaces are considered for the mixed formulation and a preconditioned conjugate gradient iterative process for solving the saddle point reformulation is proposed. We provide a preconditioning strategy that works for a large range of the perturbation parameter. Using the concepts of optimal test norm and saddle point reformulation we provide an efficient discretization strategy that works for uniform and non-uniform refinements.

ABSTRACT. Engineering problems are often characterized by significant uncertainty in their material parameters. Multilevel sampling methods are a straightforward manner to account for this uncertainty. A popular multilevel method consists of the Multilevel Monte Carlo method (MLMC). First developed by Giles, this method relies on a hierarchy of successive refined Finite Element meshes of the considered engineering problem, in order to achieve a computational speedup. Most of the samples are taken on coarse and computationally cheap meshes, while a decreasing number of samples are taken on finer and computationally expensive meshes. Classically, the mesh hierarchy is constructed by selecting a coarse mesh discretization of the problem, and recursively applying an h-refinement approach to it. This approach is referred to as h-MLMC. When comparing h-MLMC to simple Monte Carlo, the former achieves a speedup of a factor 10 with respect to the latter. However, the mesh hierarchy used for the h-MLMC method increases the number of degrees of freedom almost geometrical with increasing level, leading to a large computational cost. An efficient manner to reduce this computational cost, is by means of the novel sampling method developed by the authors, called p-refined Multilevel Quasi-Monte Carlo (p-MLQMC). The p-MLQMC method combines a hierarchy of p-refined Finite Element meshes, with a deterministic Quasi-Monte Carlo sampling rule based on rank-1 lattice rules. This combination significantly reduced the computational cost with respect to h-MLMC. However, the p-MLQMC method presents the practitioner with a challenge. This challenge consists in adequately incorporating the uncertainty, represented as a random field, in the Finite Element model. For h-MLMC, one typically uses the midpoint method to achieve this. The uncertainty is represented as a random field obtained by means of a Karhunen-Loève (KL) expansion, of which point evaluations are then computed at the element centers. These point evaluations are then assigned to the respective elements. For p-MLQMC, the integration point method is used. Here point evaluations of the random field are computed at certain carefully chosen spatial locations, and assigned to the quadrature points of the elements. These quadrature points are used to obtain the element stiffness matrices by means of numerical integration. For the p-MLQMC method, we have investigated how these evaluation points are to be selected in order to achieve a good variance reduction across the levels, and thus a lower total computational cost. By carefully selecting the evaluation points, a strong positive correlation between the solutions of two successive levels is obtained, which results in a good variance reduction across the levels. We have shown that using sets which consist of nested evaluation points across the mesh hierarchy, i.e., the Local Nested Approach (LNA), yields a speedup up to a factor 5 with respect to sets which consist of non-nested evaluation points, i.e., the Non-Nested Approach (NNA). Furthermore, we have shown that the p-MLQMC method combined with the LNA approach (p-MLQMC-LNA) yields a speedup up to a factor 70 with respected to h-MLMC. Currently, our research focus lies on implementing higher order Quasi-Monte Carlo rules. This research path shows promising results for further computational savings in the p-MLQMC method. All the aforementioned implementations are benchmarked on a slope stability problem, with spatially varying uncertainty in the ground. The goal of a slope stability problem is to assess the stability of natural or man made slopes. This is done by investigating the vertical displacement of the top of the slope. We have discretized the slope stability problem by means of an unstructured mesh consisting of triangular Lagrangian elements.

ABSTRACT. Two-dimensional, two-phase random heterogeneous microstructures are simulated using the planar Boolean Model, a standard random set model of stochastic geometry. The corresponding effective material parameters on the macroscale, computed from numerical homogenization assuming a Neo-Hookean material model on both scales, are related to the microscale parameters. The size of the RVE is considered wrt. the parameters of the stochastic model.