ABSTRACT. A parallel 3D two-phase incompressible flow solver is developed for predicting free-surface wave fields around ship hulls using PETSc. The solver is built on multiple PETSc DMDA objects with a structured-grid finite-volume discretization and a fractional-step pressure-velocity coupling. Free-surface evolution is tracked by a volume-fraction transport method. Immersed hull geometry is handled within the same DMDA framework and reflected consistently in the linear system assembly. The pressure system is solved via PETSc KSP with algebraic multigrid preconditioning. Preliminary results for surface-piercing body flows demonstrate stable wave formation and satisfactory conservation properties. Future work targets improved interface resolution and hull geometry fidelity.

ABSTRACT. The concept of an N-dimensional array, commonly referred to as ‘tensors’, is prevalent across nearly all scientific software to represent structured data. Its popularity is due to two factors: (1) it is very natural to reason about structured data in terms of axes and indices, and (2) it is straightforward to write very efficient code, across a variety of platforms, for performing tensor-like operations. Packages like PyTorch and JAX expose tensor operations as an embedded domain-specific language (DSL) in Python that they just-in-time (JIT) compile to high performance code.

It is obviously appealing to apply this type of abstraction to continuum mechanics simulations. However, tensors are an insufficient abstraction for structured data living on an unstructured mesh. Vectors of multiple fields or representing unknowns on different types of mesh entities cannot losslessly be represented as multi-dimensional tensors because the tensor dimensions are inconsistent. Furthermore, the data structures are distributed in parallel and may have ghost regions that need updating.

pyop3 is the abstraction that bridges this divide. It generalises the concept of multi-dimensional data to the tree-like structures seen with unstructured mesh computations. With pyop3, like PyTorch and JAX, users can express near-arbitrary operations over mesh-like data structures using familiar indexing syntax and execute them efficiently in parallel. pyop3 is (or will soon be) the execution engine underpinning Firedrake.

ABSTRACT. Modeling the wave-breaking phenomena is an important technology to understand the transfer of energy, momentum, and mass across the air–water interface in ocean-climate research. In this work, we consider a three-dimensional phase-field model based on the quasi-incompressible Navier–Stokes equation with the energy component, and the conservative Allen–Cahn equation, incorporating liquid–vapor phase change and volume variation. The set of PDEs is strongly coupled and highly nonlinear, making the computation of the wave breaking problem extremely expensive. To address this challenge, we develop a fully decoupled numerical scheme with a Helmholtz smoother for the surface-tension force term. To solve the resulting linear system of equations, a highly parallel solution strategy is employed, with a RAS preconditioned GMRES for the linear systems arising from the conservative Allen–Cahn, energy, and velocity equations, and a geometric–algebraic MG preconditioned CG for the pressure Poisson and Helmholtz systems. Numerical experiments demonstrate that the proposed strategy is efficient and scalable for a benchmark wave-breaking problem.

ABSTRACT. Efficient solution of incompressible flow and fluid–structure interaction (FSI) problems is challenging due to the saddle-point structure of the resulting linear systems and the strong coupling between velocity and pressure variables. Robust preconditioning strategies are therefore essential to ensure scalable convergence of Krylov subspace methods. In this work, we consider incompressible flow problems with moving boundaries formulated within an Arbitrary Lagrangian–Eulerian (ALE) framework. The governing equations are discretized in time using a backward finite difference scheme, while the nonlinear convective term is treated using a semi-implicit linearization strategy. The linearised system is then discretised in space using Taylor-Hood elements and implemented in FreeFEM++. The resulting algeraic equations are solved using PETSc. We investigate a modified Augmented Lagrangian (mAL) approach as a preconditioner within a field-split framework. The coupled velocity–pressure system is decomposed into separate blocks, where the velocity subsystem is solved using GMRES with algebraic multigrid preconditioning, and the pressure subsystem is treated using a lightweight iterative scheme. The global system is solved using flexible GMRES. The performance of the proposed mAL preconditioner is assessed on a three-dimensional broken dam problem with moving boundaries. Mesh motion in the ALE formulation is obtained by solving a Poisson equation using GMRES with hypre preconditioning. Numerical results demonstrate mesh-independent convergence and robustness with respect to mesh refinement. Parallel scalability is also investigated. To validate accuracy, we consider the classical three-dimensional lid-driven cavity flow benchmark. The results confirm that the proposed solver achieves both accuracy and efficiency for incompressible flow problems. These findings demonstrate that the mAL-based preconditioning strategy, combined with PETSc’s field-split framework, provides an effective and scalable approach for moving-boundary flow simulations, with planned extensions to fully coupled three-dimensional FSI problems.

ABSTRACT. We present a runtime-data-driven enhancement preconditioner for the incompressible Navier-Stokes equations discretized with a high-order finite element method. The method is designed for the situation that a given preconditioner is available but the resulting Krylov convergence is too slow. In such a case, vectors generated during the iteration, such as residual vectors and intermediate solution vectors, are collected and analyzed using an unsupervised learning technique. The analysis identifies the slow-convergence features that we use to construct some enhancement preconditioners that are additively incorporated into the given preconditioner. The method is implemented in the Firedrake/PETSc framework. Numerical experiments show that the resulting enhancement can reduce the number of Krylov iterations in the subsequent inexact Newton solves considerably.

ABSTRACT. petsc4py, the Python bindings to PETSc, has been tremendously successful at making PETSc more widely available to the scientific community, and without it Firedrake could not exist in its current form. However, by design, petsc4py is a faithful reproduction of the PETSc C API in Python, and so does not take advantage of any 'Pythonic' design patterns or language features to simplify the user interface.

petsctools is a suite of PETSc-related routines initially developed in Firedrake that have been pulled out into a standalone Python library (https://www.firedrakeproject.org/petsctools/). In particular it is the new home for the “OptionsManager” class that enables PETSc options to be specified from Python dictionaries as well as from the command line. Alongside this, petsctools also provides routines for examining the PETSc configuration, setting up Python preconditioners, registering citations, writing Cython code, and an “AppContext” for passing arbitrary Python types to PETSc objects, with more planned.

We believe that petsctools will be useful to anyone building Python applications on top of PETSc.

ABSTRACT. We investigate the coupling of nonlinear water-wave motion to heaving buoy wave-energy dynamics [1] in the presence of an inequality constraint. Building on augmented Lagrangian variational principles (VPs) developed by Burman [2] and others, we impose constraints of the form G(q)≥0, where q involves system variables, through a Lagrange multiplier λ. The strict Kuhn–Karush–Tucker (KKT) conditions {λ G=0,G(q)≥0, λ≤0} are replaced by smooth approximations of the involved function F(c G(q)−λ)=max(c G(q)−λ,0), with smoothing parameter c>0, allowing explicit computation of the multiplier λ as (part of) a force. The approach combines: (a) an Average Vector Field (AVF) energy-conserving time-stepping method, extended to water-wave systems with an auxiliary field, enforcing energy conservation in the discrete system; and, (b) (novel) smooth relations λ(G) that regularise the KKT conditions. This framework has been implemented and tested in Firedrake, leading to improved and surviving benchmarks for a point particle under gravity bouncing off a rigid table and forced nonlinear water waves in a horizontal channel causing buoy motion in a wave-enhancing contraction. The extended automation within Firedrake and especially subtleties in dealing with solver singularities in the AVF approach will be reported.

[1] O. Bokhove, A. Kalogirou and W. Zweers (2019) From Bore–Soliton–Splash to a new wave-to-wire wave-energy model. Water Waves 1, 217-258. [2] E. Burman, P. Hansbo and M.G. Larson (2023) The augmented Lagrangian method as a framework for stabilised methods in computational mechanics. Archives of Computational Methods in Eng. 30, 2579–2604.

ABSTRACT. Nonlinear Eigenvalue Problems with Eigenvector Dependency (NEPv) are critical in computational physics and chemistry, forming the mathematical core of simulations that calculate electronic structures and stable equilibrium states in complex systems, among many other scientific applications.

We present a new SLEPc module to solve the NEPv problem, enabling users to easily and efficiently solve this type of problem through algorithms such as the Self-Consistent Field (SCF) iteration with different types of acceleration, and leveraging SLEPc's linear eigensolvers.

ABSTRACT. Traditional numerical solvers are typically optimized for specific hardware, making them difficult to adapt to modern architectures - for example, porting legacy code to GPUs is time-consuming and often requires extensive rewriting and re-optimization. By combining domain-specific languages (DSLs) with automated code generation, the level of abstraction is raised, enabling the generation of high-performance code across diverse hardware architectures. Moreover, providing users with a high-level problem specification facilitates the development of complex PDE solvers in a form closer to continuous mathematics, reducing code complexity and maximizing reuse.

Devito, a DSL and compiler for finite-difference (FD) solvers, has been extended to integrate iterative solver functionality through an interface with PETSc, enabling application to a range of new scientific problems, such as computational fluid dynamics (CFD). As an industry-standard framework, Devito automates the generation of highly optimized explicit FD kernels and stencil computations and has been extensively used in large-scale seismic inversion and medical imaging applications. The new developments introduce automatic generation of matrix-free routines in Devito, allowing interaction with PETSc’s suite of solvers. Key enhancements include support for iterative solvers, implicit time-stepping, coupled solvers, and matrix-free preconditioning. These features are fully integrated into Devito’s symbolic API while maintaining compatibility with staggered grids, subdomains, and custom stencils. By introducing new compiler passes into Devito, we enhance its capabilities, enabling it to address a broader range of high-performance computing challenges, including incompressible flow problems in CFD.