ABSTRACT. GRS develops the system code package AC², which is used to simulate the behavior of nuclear reactors from normal operation up to design extension conditions and severe accidents with radioactive releases to the environment. PETSc is introduced to AC² via the GRS-developed library NuT (Numerical Toolkit), which provides the AC² codes with centralized and easy access to numerical algorithms and data structures based on PETSc counterparts. A continuously growing number of numerical tasks in AC² is handled by NuT. Therefore, a smooth and robust integration of PETSc is crucial to the operability of AC².

In this talk we present our approach to easily test the effects of, say, a new PETSc version on any downstream project (NuT, simulation codes). For this, we combine GitLab and CMake techniques supported by the publicly available script library CMakeIt. The approach makes use of the existing CI test environment of AC² and takes multiple platforms (MS Windows or Linux-based) and configurations into account. Furthermore, it features the amiable property not to taint the CI pipeline history of any downstream project.

ABSTRACT. PETS provides interfaces for Fortran and Python. Several Julia interfaces existed, but all only have limited functionality. Recently, we have extended the PETSc.jl package such that it now covers (nearly) the full PETSc library. There are a number of advantages compared to other attempts: - It can be very easily installed by typing "add PETSc" in the Julia package manager. - We also developed precompiled PETSc binaries and MPI for the most important operating systems (using BinaryBuilder). This is also compiled versus MPITrampoline, which allows deploying the binaries onto different HOPC systems by compiling a small wrapper library versus the local MPI installation. - We provide both a high-level and a low-level interface. The low-level interface automatically wraps nearly the full PETSc library with over 3000 functions, whereas the high-level interface gives a more Julia-like experience but is currently limited to part of the library. - Automatic testing and CI/CD are performed with currently >50’000 tests. - Users can combine features from Julia, the Julia ecosystem, such as automatic differentiation and plotting with PETSc solvers. - The resulting codes are much more compact than their counterparts in lower-level languages; yet, users still have access to all PETSc features, such as multigrid solvers for DMDA or DMStag grids. - It allows running code on both a local workstation and on a large HPC system.

Here, we will summarise some of the work done to achieve this and show scalability results of some example codes.

ABSTRACT. We consider the resolution of saddle-point linear systems arising from the discretisation of coupled or constrained systems. In many cases, such systems are challenging to solve by means of iterative methods such as GMRES and the development of efficient preconditioners is an active field of research~\cite{begoli05}. In this work, we present a robust block-preconditioner for a $2\times 2$ block-system of the form \begin{equation} \mathcal{A}= \begin{pmatrix} A&B^T\\ B&-C \end{pmatrix}\,, \end{equation} with $A$ a symmetric positive definite (SPD) block, and $C$ a positive semi-definite block.

Such problems have been extensively studied especially in academic communities and different families of preconditioners have been proposed~\cite{bewa08,Elman,Notay2017}. In the context of our study, the matrices came from the numerical approximation of the incompressible Navier-Stokes equations on general meshes for nuclear safety analyses using the finite volume method PolyMAC~\cite{IJNMFBacq,bagend23}. \\

In this talk, we highlight an issue that has not been much investigated to the extent of our knowledge. In some physical applications, the theoretical hypotheses required by the state-of-the-art linear solvers are often not satisfied. As a consequence, even the most simple systems with a matrix in the form of Equation~(1) can become challenging to solve. This is the case for instance when the diagonal dominance of the block $A$ is lost. Indeed, one typical classical approach is to use a block diagonal preconditioner for the matrix~(1) \begin{equation} P_{Class}= \begin{pmatrix} D_A&0\\ 0&-C-BD_A^{-1}B^T \end{pmatrix}\,, \end{equation} where $D_A$ the diagonal part of $A$ is used as an approximation of $A$, and $-C-BD_A^{-1}B^T$ is used as an approximation of the schur complement $-C-BA^{-1}B^T$. However for some numerical methods, the lack of diagonal dominance of $A$ result in~(2) yielding an inefficient preconditioner.\\

Such a difficulty occurs for instance on distorted meshes. Industrial solvers then resort to direct solvers with all the ensuing limitations. To tackle this difficulty, we present an algebraic preconditioner with increased robustness for systems such as Equation~(1). The key idea involves an algebraic transformation of the system (see \cite{Notay2017,IJNMFBacq}) that compensates the loss of diagonal dominance : \begin{equation} \hat{\mathcal{A}}= \begin{pmatrix} A&B^T\\ B&-C \end{pmatrix} \begin{pmatrix} \mathbb{I}&-D_A^{-1}B^T\\ 0&\mathbb{I} \end{pmatrix} = \begin{pmatrix} A&-(\mathbb{I}+AD_A^{-1})B^T\\ B&-C-BD_A^{-1}B^T \end{pmatrix} \,. \end{equation} The new linear system with matrix~(3) can then be solved efficiently using the following block diagonal preconditioner \begin{equation} P_{Trans}= \begin{pmatrix} 2D_A&0\\ 0&amg(-C-BD_A^{-1}B^T) \end{pmatrix}\,, \end{equation} where $amg(-C-BD_A^{-1}B^T)$ denotes any classical algebraic multigrid approximation of $-C-BD_A^{-1}B^T$ (note that $C+BD_A^{-1}B^T$ is SPD).\\

Since the approach we propose is algebraic in nature, it can be applied to a broad class of problems beyond nuclear safety, although we mostly tested it on problems from fluid dynamics. We implemented and tested several variants (left transform, right transform, symmetric transform) in C using PETSc 3.25 \cite{PETSc}. The MATNEST structure enabled us to manipulate each block to build the transformed matrix~(3) at a very low cost. The PCFIELDSPLIT structure enabled us to build block preconditioners, and the PCGAMG structure provided us with an efficient multigrid preconditioners. The resulting code is open source, available online \cite{SddlePointGithub} and some illustrative scripts are meant for eventual integration in PETSc tutorial examples.\\

We first highlight in the talk the loss of robustness of classical approaches when considering distorted meshes. In a second time, we describe our preconditioning approach based on an algebraic transformation of the matrix~(1). In the third part we give some numerical results that show the impressive convergence of the new preconditioner on problems of industrial complexity. As for the perspectives we propose an international benchmark on saddle point problems with meshes and matrices originating from fluid dynamics simulations \cite{Benchmark}.

\begin{thebibliography}{9.} \frenchspacing \bibitem{begoli05} M. Benzi, G. H. Golub, J. Liesen. Numerical Solution of Saddle Point Problems. Acta Numerica, 14, 1–137, 2005. doi :10.1017/S0962492904000212. \bibitem{bewa08} M. Benzi, A. J. Wathen. Some Preconditioning Techniques for Saddle Point Problems, pp. 195–211. Springer Berlin Heidelberg, Berlin, Heidelberg, 2008. doi :10.1007/978-3-540-78841-6\_10. \bibitem{Elman} Howard C. Elman, David J. Silvester, and Andrew J. Wathen. Finite elements and fast iterative solvers: with applications in incompressible fluid dynamics. Oxford University Press, 2014 \bibitem{Notay2017} Y. Notay, Algebraic Multigrid for Stokes Equations, SIAM J. Sci. Comput., 39 (2017), pp. S88–S111. \bibitem{IJNMFBacq} Bacq, P.-L., Gerschenfeld, A. and Ndjinga, M. (2025), PolyMAC: A Staggered Finite Volume Method on General Meshes for Incompressible Navier-Stokes Equations. Int J Numer Meth Fluids. https://doi.org/10.1002/fld.5398 \bibitem{bagend23} P.-L. Bacq, A. Gerschenfeld, M. Ndjinga. PolyMAC : Staggered Finite Volume Methods on General Meshes for Incompressible Navier–Stokes Problems. In E. Franck, J. Fuhrmann, V. Michel-Dansac, L. Navoret, eds., Finite Volumes for Complex Applications X—Volume 1, Elliptic and Parabolic Problems, pp. 149–156. Springer Nature Switzerland, Cham, 2023. \bibitem{PETSc} Satish Balay et al., {PETS}c {W}eb page, \url{https://petsc.org/}, 2026 \bibitem{SddlePointGithub} M. Ndjinga, Saddle point linear solver repository, \url{https://github.com/ndjinga/SaddlePointLinearSolver} \bibitem{Benchmark} M. Ndjinga, NavierStokes 3D Benchmark repository, \url{https://github.com/ndjinga/NavierStokes3DBenchmark}, 2023 \end{thebibliography}

ABSTRACT. Typical domain decomposition experiments use as many MPI processes, and typically also as many CPU processor cores, as subdomains, each of the CPU cores corresponding to one MPI process and handling one subdomain. In the context of two-level Schwarz domain decomposition methods, and using the -pc_asm_blocks PETSc option, we here explore the potential benefits of using more subdomains than CPU cores, i.e., having each core, corresponding to one MPI process, handle more than one subdomain, and we show how this can be beneficial. We provide results in terms of time-to-solution, memory usage, and energy-to-solution, the latter metric being of increasing importance in supercomputing centers. Moreover, we also explore the possibilities offered in this context by the recently introduced -pc_precision PETSc option, having the local and/or coarse solves handled in single (instead of double) precision by MUMPS.

ABSTRACT. We present recent mathematical and software developments that enable robust and scalable solvers in Firedrake through PETSc. We apply the theory of Schöberl on multigrid for nearly-indefinite systems to both conforming macroelements and non-conforming elements. Two crucial components of the multigrid solver are the relaxation and the inter-grid transfer operators, both of which need a kernel-capturing property.

For the relaxation we employ the additive Schwarz method on small patches constructed around topological mesh entities. This ensures a robust relaxation with respect to the problem parameters and the discretization. The key development is the grouping of the patches into non-overlapping subsets via a coloring of the mesh entities. While producing a mathematically equivalent space decomposition, the mesh coloring significantly improves performance by reducing overheads in the setup and application of the individual patch solvers.

For the transfer operator we implement two approaches based on standard interpolation into the skeleton of the coarse or the skeleton of fine mesh, followed by harmonic extension into the interior. We compare the two approaches: the one on the fine skeleton requires a more expensive interpolation using a composite quadrature rule to evalute degrees of freedom, whereas the one on the coarse skeleton involves a more expensive harmonic extension.

We demonstate our approach with numerical examples from nearly incompressible elasticity. We also apply these solvers to adaptively-refined mesh hierarchies. These hierarchies are constructed with Netgen and ngsPETSc by estimating the error of a user-specified goal functional.

ABSTRACT. Domain decomposition methods are a cornerstone of scalable solvers for large sparse linear systems, and HPDDM provides efficient implementations exposed in PETSc through PCHPDDM. Previously, HPDDM relied on its own backend, limiting its ability to leverage GPU acceleration available in PETSc. In this work, the HPDDM backend is migrated to PETSc-native data structures and operations, enabling seamless GPU support and improved performance portability. I will also highlight additional contributions to both HPDDM and PETSc that improve interoperability, code structure, and extensibility.

ABSTRACT. PFLARE is an open source library available through the PETSc configure which aims to provide fast & scalable iterative methods for asymmetric linear systems, in parallel and on both CPUs and GPUs.

This talk focuses on solving advection equations, as traditionally the only scalable approach in parallel involves using sweeps (i.e., matrix-free Gauss-Seidel methods) on structured grids. Multigrid methods are now available that target asymmetric linear systems, specifically those which construct approximate ideal restrictors (AIR).

We present results from using an AIR variant known as AIRG, which uses GMRES polynomials to form approximate inverses and has been designed for use in parallel and on GPUs. We use AIRG to solve advection equations on Lumi-G, a 380 PFLOP GPU machine and demonstrate good weak scaling on both structured and unstructured grids.