ABSTRACT. A range of strategies have been developed which reduce the storage requirements of reverse mode ("adjoint mode") algorithmic differentiation. These trade storage for computation: replacing complete storage of all variables needed by the adjoint with some combination of storage and recalculation. The implementation of efficient checkpointing strategies is one of the major challenges in applying algorithmic differentiation in practice.

This talk will detail the challenges which are encountered when applying algorithmic differentiation to tools which apply a 'dynamic computational graph' approach (the approach used by pyadjoint). Generalization to high order adjoint calculations will be discussed.

Recently new checkpointing approaches have been developed which have improved performance -- and which for example out-perform the well-known revolve algorithm. This talk will outline the implementation of one such scheme, demonstrating that improved performance is achievable, in practice, with codes written using Firedrake.

ABSTRACT. The use of adjoint calculations to compute the gradient of a quantity of interest resulting from solving a system of partial differential equations (PDEs) is widespread and well-established. In Firedrake, the adjoint-based gradient is computed using algorithmic differentiation and has been successfully applied in several scenarios, such as optimisation problems. For such cases, it is common to define a reduced functional objective, which enables the computation of an adjoint-based gradient and the subsequent recomputation of the functional. Firedrake has introduced the Ensemble Reduced Functional, a novel object that allows the simultaneous computation of functionals and their gradients. A direct application of the Ensemble Reduced Functional is in the seismic inverse problem, where multiple functionals are associated with different wave sources. In this case, the total functional is the summation of the wave sources' functionals. The Ensemble Reduced Functional enables simultaneous execution of these wave sources while allowing for mesh parallelisation, resulting in faster simulations. However, this approach requires high memory usage, as adjoint PDEs are executed in reverse time and depend on forward data. Storing the entire forward state for the adjoint calculation results in a memory footprint that grows linearly with the number of time steps, potentially exhausting available computational resources. Checkpoint algorithms enable the management of memory usage in Firedrake. These algorithms store only the state required to restart the forward calculation from a limited set of steps, effectively controlling the memory footprint. In summary, this presentation aims to demonstrate the usage of these new features available in Firedrake for executing gradients, hopefully enabling faster simulations while managing memory more effectively in various applications.

ABSTRACT. Numerical models are critical tools in understanding the behaviour of coastal ocean systems. In particular, tidal stream energy deployments are set to expand in capacity by over 1200% in the next 4-5 years and thus determining the spatial distribution of the resource is warranted. An accurate calibration of coastal models can reduce uncertainty and improve confidence for investors. In shallow-water equation models, the field typically used to improve agreement with measured values is the bottom friction, which invariably must account for broader physical and numerical dissipation processes that arise at unresolved spatiotemporal scales. This has ordinarily been performed through enumeration of a set of reasonable parameters, assuming a uniform value across the domain. However, both skin friction and form drag will vary spatially and therefore other representations have been proposed including simple spatial segmentation, mapping based on bed sediment particle size data, and an 'independent points scheme'. Resource models are typically only validated with simple sensitivity approaches to provide high-level assessments, but wider-use coastal ocean models have been formally calibrated based on altimetry or tidal gauge (depth) data using adjoint-based techniques, Bayesian inversion and Kalman filtering. However, tidal stream deployments are sensitive to the velocity field which is much more challenging to reproduce accurately. Thus, numerical models now need to be calibrated based on Acoustic Doppler Current profiler data to allow layouts to be correctly micro-sited. This requires the development of effective strategies to frame the optimisation problem. This talk investigates the use of adjoint-based optimisation of the bed friction field in the Thetis coastal ocean model considering various forms of friction representation. This begins with idealised set-ups using the various methods of friction representation, increasing in complexity to the site of an upcoming tidal array deployment. This talk will focus on the progress made and the challenges faced in using this method.

ABSTRACT. In this study, we develop a digital twin of a Westerly granite rock sample using pyadjoint in the Firedrake framework to investigate the evolution of fluid transport properties and heterogeneity during fault growth. Our goal is to provide detailed insights into the dynamics of fluid flow during fault development through advanced numerical simulations informed by experimental data.

We conducted a series of experiments where a Westerly granite sample was subjected to triaxial conditions, inducing quasi-static failure. Fault growth was closely monitored using acoustic emission data, and periodic flow-through tests were performed to measure pore pressure heterogeneity and flow rates. Using these experimental observations, we numerically reconstructed the exact geometry of the tested sample and employed three-dimensional finite element simulations of Darcy flow within an adjoint framework. This allowed us to estimate the heterogeneous fluid flow properties at various stages of the faulting process through least-square minimisation.

Our digital twin enabled us to identify critical stages in the faulting process where permeability undergoes significant changes. Notably, following peak stress, the permeability of the fault zone increased up to 150 times that of the bulk material. A subsequent significant increase, up to 400 times the bulk permeability, was observed when the fault slip ranged between 0.6 and 0.7 millimetres. Beyond these stages, no substantial permeability changes were detected. Additionally, we quantified the permeability heterogeneity along the shear fault zone, identifying variations up to 8 times between different fault regions, depending on the specific faulting stage.

This work underscores the potential of digital twins in advancing state-of-the-art rock mechanics experiments, offering valuable insights into the understanding of geophysical properties like the evolution of fluid transport during fault growth.

ABSTRACT. Recently, there has been a lot of interest in machine learning based surrogate models for atmospheric fluid dynamics. Successful approaches such as GraphCast [Lam et al, Science, 382(6677), 2023] use an encoder-processor-decoder architecture where the processor describes the propagation of information on a Graph Neural Network (GNN) via message passing. Since the dynamics are approximated on a low-dimensional space, these models are much faster than traditional methods that solve the underlying PDE numerically on a fine mesh.

We explore a variation of this architecture where the processor is replaced by the solution of a time-dependent differential equation in a low-dimensional latent space. Compared to GNNs, this allows the application of standard techniques from numerical analysis, thus potentially improving stability and interpretability for example by exactly enforcing conservation laws.

To implement this, we use the recent Firedrake/PyTorch interface [Bouziani & Ham, 2023] to train and solve time-dependent PDE surrogate models on the sphere. The encoder in our model combines the interpolation of the initial condition to the latent space on a vertex-only-mesh with a learnable embedding; the decoder has a similar structure based on the adjoint of the interpolation. Our model is trained on numerical solutions of PDE which have been calculated a-priori using Firedrake. The aim of our work is to combine the reliability of Finite Elements with the efficiency of Neural Networks surrogates to produce a competitive model which has the potential to be applied to time-critical applications in weather forecasting.

ABSTRACT. Neural operators have been very effective in learning solution operators associated with partial differential equations. However, physical systems with complex geometries or boundary conditions or the discretization choice of the input-output functions are often challenging tasks for existing approaches. We introduce a structure-preserving operator learning framework that leverages structured finite element spaces to learn operators associated with partial differential equations. Our approach is easy to implement using Firedrake, has theoretical guarantees and automated structure-preserving spatial discretizations, is end-to-end differentiable, and preserves the mathematical and physical properties of the PDE.

ABSTRACT. After a mini introduction to shape optimization and firedrake's shape optimization toolbox fireshape, I will present fireshape's latest developments and discuss future plans.

ABSTRACT. A magnetic cloak is a device that hides objects inside and makes them undetectable to magnetic fields outside. Current cloaking designs face two main challenges:(1) a slow development and testing process, and (2) applicability limited to simple geometries.

To overcome these limitations, we propose a fast and accurate numerical approach for achieving magnetic cloaking on arbitrary geometries. This approach involves solving a curl-curl PDE-constrained optimization problem using Firedrake with pyadjoint.

ABSTRACT. We present our finite element solutions to a continuum model of vascular tumour growth. The model describes the evolution of a multiphase system of four incompressible phases – healthy cells, tumour cells, blood vessels and extracellular material – via three systems of partial differential equations (PDEs). These systems comprise: mass balance equations for the volume fraction of each phase; momentum balance equations for the velocity and pressure fields of each phase; and a reaction-diffusion system describing the transport and consumption of two essential nutrients, oxygen and glucose, that feed the growth of the tumour phase.

We employ a carefully-selected combination of continuous and discontinuous finite element spaces to discretize the 3D spatial domain, coupled with the three-stage strong-stability-preserving Runge-Kutta scheme for the time stepping. Well-designed preconditioners are used to accelerate the linear algebra required at each time step significantly. A range of numerical results will be presented for discussion. Our implementation is in python, making use of the Firedrake framework.

ABSTRACT. Silicon is a frequently used active material in the negative electrode of lithium ion batteries as it provides substantial improvements in the energy density. Due to large volume changes during cycling, the Si content in state-of-the-art electrodes is typically rather low. As significantly higher Si contents are desirable the effects of structural changes and electrolyte displacements have to be analysed. Accurate models will enable the study of possible local electrolyte depletions in order to mitigate degradation processes and improve cycle life of the batteries.

In our work we aim to extend the mesoscopic transport theory of Li-ion batteries by Latz et al. [1] including single phase flow through the porous electrode media. This Darcy flow is generated by the change in Si volume during battery operation. After eliminating the velocity, a discontinuous Galerkin (piecewise constant) formulation equivalent to the finite volume method with upwind schemes is obtained [2]. First toy problems were implemented in Firedrake and tested for convergence properties. Next step is the realisation of the full cell geometry incorporating electrochemical transport. Future work includes the consideration of a radial pseudo dimension within each electrode discretization voxel to account for the diffusion of Li within individual spherical particles and numerical optimization e.g. through stabilizing terms or adaptive mesh refinements.

[1] Arnulf Latz et al 2011, Journal of Power Sources 196 3296-3302 [2] Thomas Roy et al 2019, Journal of Computational Physics 395 636–652

ABSTRACT. The transformation of the European energy sector to more efficient and sustainable technologies is a matter of utmost importance in the face of climate change. However, considering the fact that renewable energy sources are inherently volatile, the development of novel and innovative electrochemical storage systems “calming the waves” of volatility has become a very active research area. Along with these efforts and the prospect of large electric energy storages, the interest in battery-powered aerospace applications has also grown. The full potential of novel battery types in one of the most demanding engineering environments, in which energy, power, and weight requirements must be jointly met at the same time, is yet still unclear.

Lithium-sulfur batteries as conversion type batteries are currently believed to be one of the most promising candidates due to their exceptionally high theoretical gravimetric energy density, being just a quarter of kerosene burnt in civil aero-engines. However, there is already empirical evidence that building cells with high gravimetric energy density and power density is conflicting due to kinetic limitations and, hence, contrary to kerosene-based propulsion. To exactly pin down this optimum, which is also believed to depend on the microstructure of the porous cathode, is challenging, even for cutting-edge experimental operando methods. Therefore, we have developed a 3D scale-resolving Firedrake framework as a complementary “numerical operando technique”, which aims at narrowing the optimization window and helping to find the most relevant parameters controlling the multi-scale transport process in the porous structure. The locally coarse-grained continuum model is currently solved using a DG0 approach in combination with an adaptive time-stepping procedure, which will be presented alongside with first vivid insights into the nonlinear transport process.

ABSTRACT. This talk will introduce a novel construction of adaptive and yet structure preserving finite element discretizations with function spaces induced by subdivision. In many applications, for example in geophysical fluid dynamics, adaptive and structure-preserving methods can be highly beneficial to simulate the long-term evolution of a multi-scale system with several invariants of motion like the total energy. If the discretizations of such systems do not preserve these invariants, the simulation results can differ significantly from the true physical behaviour of the systems. Combining the benefits of structure preservation and adaptive finite elements is notoriously difficult. If no special care is taken, adaptive mesh refinement algorithms of standard finite element approaches usually lose the property of structure preservation. Alternatively, the refinement can be chosen to be conforming, which in turn leads to unnecessary propagation of the refinement because surrounding cells need to be refined as well. On the other hand, IGA tensor-product techniques suffer from mesh topology restrictions. For this reason, we chose to build our function spaces upon subdivision. We extend the work of [1], who introduced vector field subdivision schemes that commute with the standard vector calculus operators like the gradient or the curl. Translating their work to the finite elements realm yields de-Rham-complex-preserving finite elements for scalar functions, vector fields, and density functions. We added adaptivity to their structure-preserving discretization by leveraging the hierarchy of the basis functions induced by the subdivision algorithm. By carefully keeping track of the introduced degrees of freedom across the refinement levels, we maintain a discrete de Rham complex and thus enable structure-preserving simulations. Our method was verified by simulating the Maxwell eigenvalue problem, a well-known test case that reproduces the analytical eigenvalues if the chosen finite element spaces constitute a discrete de Rham complex. We show that our discretization indeed yields the correct spectrum and investigate the computational effort and accuracy gains of our method.

[1] Ke Wang, Weiwei, Yiying Tong, Mathieu Desbrun and Peter Schröder: “Edge subdivision schemes and the construction of smooth vector fields”. In: ACM Transactions on Graphics (2006)

ABSTRACT. We specialize in developing high resolution Digital Light Processing Stereolithographic (DLP-SLA) 3D printing for microfluidic (lab-on-a-chip) device fabrication. Fabrication of such devices focuses on negative rather than positive features because in applications fluid to be processed and manipulated in the device will occupy the negative features. To achieve the necessary negative feature resolution, we have built our own custom 3D printers and resins, and have developed all of the required control software and fabrication processes. We routinely fabricate microfluidic devices that integrate both passive and active (valves, pumps, and mixers) components on a size scale that so far has not been matched. As an example, we can now 3D print channels with cross sections as small as 1.9 microns x 2.0 microns.

In our quest for ever higher resolution, we are at the point where modeling the photopolymerization process would be useful to elucidate ultimate resolution limits and optimization strategies. To that end we have implemented a finite element approach with Firedrake to solve a set of four coupled partial differential equations (PDEs). One of the equations is directly nonlinear in the field quantity (reactive species concentration) and two of them have nonlinear reaction rates. In addition, the diffusivity in the diffusion terms in all four PDEs are nonlinear. In this presentation I will cover our solution approach, results so far, and introduce a layer-by-layer method to sequentially simulate the photopolymerization of one layer, transfer the result to a new mesh initialized with a new layer of resin, and then simulating the photopolymerization of the combined (stacked) layers, which can be extended to simulating arbitrary numbers of layers. I will also discuss pitfalls and problem areas.

ABSTRACT. Metastable alloys are known for their optimal combination of strength and ductility, attributes that are achieved through complex thermomechanical treatments involving multiple phase transformations. A particularly significant transformation is the beta-to-omega phase transformation, during which nano-sized omega phase particles are formed. To better understand this process, we developed a phase-field model capable of simulating the formation and evolution of the omega phase. Our model accounts for the transformation of the beta phase into four distinct variants of the omega phase, resulting in a system with seven global variables (three displacements and four order parameters). Achieving high-resolution simulations of omega particles within the beta matrix requires a fine computational mesh. This model is implemented using the Firedrake finite element package, enabling the use of a geometric multigrid method to solve this large-scale problem efficiently.

ABSTRACT. This presentation will explore an advanced particle filtering framework that combines tempering, jittering, and global nudging. The entire filtering process is executed with the power of ensemble parallelism, ensuring efficient and scalable computations. The effectiveness of this framework is demonstrated through its application to the stochastic Kuramoto-Sivashinsky (KS) and Euler models, showcasing its capability to handle complex, nonlinear dynamical systems.

ABSTRACT. We present contributions toward a finite element method for elliptic problems, whether the scalar Poisson problem or the elasticity problem, that incorporates experimental or synthetic constitutive data without relying on a smooth constitutive relation-between fluxes and gradients in the former case, or to stresses and deformations in the latter. This approach is known as the Data-Driven Computational Mechanics paradigm, which can be seen as a form of unsupervised machine learning applied to continuum mechanics. The ultimate goal of this work is to formulate a discrete-continuous optimization problem that seeks, under certain constraints, to assign to each material point a point in the space of fluxes-gradients (or stresses-deformations). This assignment minimizes the distance induced by some norm between the dataset and the subspace of corresponding compatible fields in equilibrium. Given the challenging nature of this optimization problem, we first devise strategies to initialize the data distribution over the domain and then solve the corresponding optimality conditions. In this presentation, we show progress in this direction, provide details of the Firedrake implementation, and illustrate with several 2D numerical examples.

ABSTRACT. We investigate flow in descending aorta with Navier's slip boundary condition. Given some measured flow data, for example modern 4D-PC MRI image, an variational data assimilation approach is used to estimate the Navier's slip parameter on the wall and the inflow velocity profile. This is achieved by using Firedrake and pyadjoint frameworks and tested on artificially generated data. We investigate some aspects of the robustness with respect to discretization and the noise level in the data.