ABSTRACT. Mixed-dimensional partial differential equations (PDEs) are equations coupling unknown fields defined over domains of differing topological dimension. Such mixed-dimensional PDEs naturally arise in a wide range of fields including geology, bio-medicine, and fracture mechanics.

Finite element discretizations of mixed-dimensional PDEs involve nested meshes of heterogeneous topological dimension. The assembly of such systems is non-standard and non-trivial especially with regard to the terms involved in the interactions between the different domains. In other words, automated solution of mixed-dimensional PDEs requires the design of both generic high level software abstractions and lower level algorithms.

In this talk, we present a set of features developed as part of an automated framework dedicated to such problems in the FEniCS Project framework. The introduced tools will be illustrated by concrete examples of applications highlighting their relevance for model in biomedicine.

ABSTRACT. Checkpointing, i.e. saving and reading results of finite element computation is crucial, especially for long-time running simulations where execution is interrupted and user would like to restart the process from last saved time step. On the other hand, visualization of results in thid-party software such as ParaView is inevitable. In the previous DOLFIN versions (2017.1.0 and older) these two functionalities were strictly separated. Results could have been saved via HDF5File interface for later computations and/or stored in a format understood by ParaView - VTK’s .pvd (File interface) or XDMF (XDMFFile interface). This led to data redundancy and error-prone workflow. The problem essentially originated from incompatibilities between both libraries, DOLFIN and ParaView (VTK). DOLFIN’s internal representation of finite element function is based on vector of values of degrees of freedom (dofs) and their ordering within cells (dofmap). VTK’s representation of a function is given by it’s values at some points in cell, while ordering and geometric position of these points is fixed and standardised within VTK specification. For nodal (iso- and super-parametric) Lagrange finite elements (Pk , dPk ) both representations coincide up to an ordering. This allows to extend XDMF specification and introduce intermediate way of storing finite element function - intrinsic to both, ParaView and DOLFIN. The necessary work was done as a part of Google Summer of Code 2017 project Develop XDMF format for visualisation and checkpointing, see https://github.com/michalhabera/gsoc-summary. New checkpointing functionality is exposed via write checkpoint() and read checkpoint() methods.

ABSTRACT. We describe ongoing work on extending the Diderot programming language for scientific visualization to handle functions produced by FEniCS. Existing FEM visualization tools limit the FEniCS user's ability to create high quality results either because they support only approximations of higher-order elements, or because they do not support writing new visualization algorithms. In contrast, the Diderot language supports higher-order elements and offers a concise and mathematically expressive language for implementing visualizations. We show how Diderot can express commonly available basic visualizations, while also allowing users to create new visualizations that accurately depict the structure of FEM solutions.

ABSTRACT. We have developed a domain-specific compiler to visualize data for the FEniCS community. In another submission we demonstrate the results of this work; visualizations of finite element data. That submission will demonstrate images that have been created and present new informed visualization techniques that can been created using the language. In this abstract submission we answer the question of “how”. We describe the development of our language and how we can support the type of high-level analysis necessary to visualize finite element data. In order to enable interesting visualizations of FEM field, it was necessary to create a compiler that understood FEM fields. Our Diderot compiler is designed to represent FEM fields internally, evaluate these fields at points, and supports a rich level of tensor operators on them. By pushing domain knowledge inside our compiler we are able to take derivatives which supports a richer class of visualization algorithms. Lastly, by developing a domain-specific compiler we can find common computations involved with field evaluation that were obscured by our previous method.

ABSTRACT. The Rapid Optimization Library (ROL) is a C++ library developed at the Sandia National Labs for large scale optimisation problems.

ROLs key features are: + Matrix-free APIs: enable direct use of Dolfin data structures and PETSc for linear algebra operations. + State-of-the-art algorithms for unconstrained optimization, constrained optimization and optimization under uncertainty; enable inexact and adaptive function evaluations and iterative linear system solves. + Special APIs for simulation-based optimization - enable a streamlined embedding into engineering applications, rigorous implementation verification and efficient use. + Modular interfaces throughout the optimization process - enable custom and user-defined algorithms, stopping criteria, hierarchies of algorithms, and selective use of a variety of tools and components.

We have written Python bindings for ROL using the pybind11 library. Furthermore pyadjoint was updated and is now able to use ROL to solve optimisation problems.

This work was partially funded by a EPSRC UK-USA travel grant.

ABSTRACT. In this talk, we will demonstrate our Julia wrapper for the FEniCS Finite Element Toolbox, part of the JuliaDiffEq suite, and its internals, along with snippets of code demonstrating various functionalities which allow the use of FEniCS in Julia. We then proceed to show how the wrapper can be extended to include further objects and functions as required. Finally, we will demonstrate some code solving the FEniCS tutorials, and show any differences in the API.

ABSTRACT. We study two techniques for correcting the geometrical error associated with domain approximation by a polygon. The first was introduced some time ago and leads to a nonsymmetric formulation for Poisson's equation. We introduce a new technique that yields a symmetric formulation and has similar performance. We compare both methods on a simple test problem.

ABSTRACT. We discuss the Firedrake implementation of a discontinuous Galerkin finite element method for the approximation of solutions to elliptic partial differential equations in nondivergence form. The domains under consideration are assumed to be piecewise sufficiently smooth and Lipschitz continuous, with piecewise curved boundaries. We consider both Dirichlet and oblique boundary-value problems.

ABSTRACT. The classical Biot and Terzaghi soil models, describing flow through a single fluid network in a porous and elastic medium, were generalized to equations describing multiple fluid network poroelasticity (MPET) by Barenblatt and Aifantis. The MPET equations have been utilized in geomechanics to simulate multiply-fractured strata for a few decades, but are now also beginning to find their application in biomechanics. Indeed, the multiple network poroelasticity theory aptly models the multiple fluid networks encountered in e.g. the brain: such as extracellular spaces, vasculature and paravasculature.

For the MPET model, the characterizing parameters for each fluid network are the Biot-Willis coefficient, the storage coefficient, and the hydraulic conductivity tensor. Each network equation also contains transfer terms which characterize the exchange in/out of the network. In the context of the brain, the storage coefficients and permeability tensors may typically be very small and with possible jumps. These regimes pose challenges for standard finite element discretizations of the MPET equations. This fact has spurred recent interest in developing numerical methods, and preconditioners, which are robust in the intrinsic limits of elastic, and poroelastic models.

A Stokes-Biot stability condition for Biot's consolidation model has shown recent computational promise of robustness in the vanishing permeability limit. In this talk we present an H(div)-based conformal mixed method for the MPET equations in the framework of Stokes-Biot stability. The method is based on a 3-field formulation via a fluid flux for each permeable fluid network. Early results concerning robustness properties of the approach in the low permeability and hydraulic conductivity limits will be discussed.

ABSTRACT. In an attempt to reconcile the advantages of a Lagrangian method and an Eulerian method, hybrid particle-mesh methods make combined use of Lagrangian particles and an Eulerian background mesh. Typically, the particles are used in the kinematic part of the problem for tracking (multi-)material quantities such as mass or momentum. The Eulerian mesh on the other hand, can be efficiently used for computing the dynamic part of the problem at the continuum level. A coupling between the Lagrangian stage and the Eulerian stage is established by means of a projection from the particles to the mesh, and \textit{vice versa}. Despite many successful applications, it is not surprising that such an approach tends to have difficulties in unifying accuracy and exact conservation, as information is repeatedly projected back-and-forth between a (discrete) set of Lagrangian particles and an Eulerian mesh field. \\ % This contribution presents a particle-mesh interaction strategy on arbitrary meshes which fundamentally overcomes aforementioned issue. Key to the presented approach is to formulate the particle-mesh interaction in terms of a PDE-constrained minimization problem in such a way that, from a mesh-perspective, the transported Lagrangian particle field ($\psi_p$, say) satisfies an advection equation. Hence, the constrained particle-mesh projection can be loosely phrased as % \begin{align} \label{eq:scalar_objective} \begin{split} &\min_{\psi_h} J = \sum_{p }^{} \frac{1}{2}\left( \psi_h(\mathbf{x}_p(t), t) - \psi_p(t)\right)^2 \\ &\text{such that a scalar advection problem is satisfied in a weak sense.} \end{split} \end{align} % with $\psi_h$ the state variable to be approximated in an appropriate function space. Furthermore, the summation over $p$ runs over all particles in the domain of interest.\\ % Starting from Eq.~\eqref{eq:scalar_objective}, we will derive the optimality system for the constrained particle-mesh projection and~-~by introducing a facet-based control variable~-~we will show that a hybridized Discontinuous Galerkin (HDG) framework naturally provides the ingredients required for the optimality control. Several theoretical results of the resulting particle-mesh projection are presented, including consistency and local conservation. An efficient solution strategy using \textit{static condensation} is proposed, and we show how this can be implemented using the FEniCS-framework. \\ % The presented approach can be of value to a range of applications. In particular, we show how the scheme can be applied in a particle-mesh operator splitting scheme to the advection-diffusion equation and the incompressible Navier-Stokes equations. Various numerical examples demonstrate the model capabilities in terms of accuracy and robustness. A high Reynolds-number test gives further evidence for the robustness and low artificial dissipation of the scheme.

ABSTRACT. The primary concern of the presentation is geometric inverse problems governed by hyperbolic partial differential equations, meaning we are interested in reconstructing geometric objects such that they reproduce a measured echo of a scanning wave. On the one hand, this means we consider automatisms to include the computation of shape derivatives, i.e., differentiating with respect to the input mesh, into FEniCS.

On the other hand, we also study non-smooth problems that arise naturally when objects with kinks are to be reconstructed. To this end, we consider using Fenchel Duality and Raviart--Thomas spaces for Total Variation denoising for images on surfaces. We conclude with numerical examples where FEniCS is interfaced with 3D scanners to conduct denoising of real world objects.

ABSTRACT. In this talk, we present a recent addition to dolfin-adjoint: support for Tensorflow models. Tensorflow is a numerical computation library which is mostly used to conduct machine learning (ML) and deep neural networks research, but which can also be applied in various other domains. With Tensorflow support in dolfin-adjoint, users can now implement mixed PDE-ML models, where the PDEs are solved with FEniCS and the ML functions are computed with Tensorflow. The high-level differentiation capabilities of dolfin-adjoint automatically compute derivatives of these models, which may be used to optimise (train) model parameters from data.

The implementation works as follows. During model execution, an algorithmic differentiation (AD) tool records a computation graph. In this graph, nodes represent high-level mathematical operations such as a PDE solve using FEniCS, or a ML function evaluation using Tensorflow. Edges in this graph represent variables, such as finite element functions or tensors, that flow between operations. Further, we implemented basic transfer functions to map tensors in Tensorflow to FEniCS data structures. From this computation graph, dolfin-adjoint derives the associated adjoint graph, in which each node is replaced by its assoiated adjoint. If a graph node represents a PDE solve, the AD capabilities in FEniCS is used to obtain its adjoint. Otherwise, if a node represents a Tensorflow operation, then Tensorflow's internal AD tool is used to obtain the adjoint version of that operation.

ABSTRACT. Glioblastoma (GBM) is the most frequent malignant brain tumour in adults. Its growth is characterized by infiltration of surrounding healthy tissue, and the formation of a necrotic core. GBM presents with varying degree of mass-effect which results in healthy-tissue deformation, midline shift or herniation. Biomechanical forces, such as those resulting from displacive tumour growth, shape the tumour environment, contribute to tumour progression and may affect treatment response and outcome.

To investigate the role of tumour mass-effect for tumour evolution, we have previously developed a mechanically-coupled reaction-diffusion model that captures three dominant aspects of macroscopic GBM growth: (a) tumour cell proliferation, (b) the diffuse invasion of the growing tumour into surrounding healthy tissue, and (c) the resulting mass effect.

Here we present an implementation of this model in FENICS and first steps towards an image-based optimization approach, based on dolphin-adjoint, to estimate patient-specific parameters from clinical magnetic-resonance imaging (MRI).

ABSTRACT. In this contribution we present a framework for expressing finite element methods on arbitrarily many intersecting meshes: multimesh finite element methods. These methods enables the use of separate meshes to discretize different parts of a computational domain that are naturally separate, such as the components of an engine, the domains of a multiphysics problem, or solid bodies interacting under the influence of forces from surrounding fluids or other physical fields. Multimesh finite element methods are also suited to problems in which the computational domain undergoes large deformations as a result of the relative motion of the separate components of a multi-body system. Such problems usually requires mesh regeneration which can be a costly and difficult procedure. We will also briefly mention implementation aspects in FEniCS.

ABSTRACT. In this presentation, we will introduce a cut finite element method for two-phase flows. In our scheme, which is based on CutFEM [1,2] the interfaces between two liquids are described by a level set function whose zero level set cuts arbitrarily through a fixed background mesh.

In intersected elements, the physical boundaries of each liquid phase are respected by integrating the weak formulation only over parts of the elements that are covered by the respective fluids. Additionally, intersected elements are doubled to allow for discontinuities in the pressure and for kinks in the velocity solution inside elements. Our method is stable independent of the cut location due to so-called ghost penalties [3], which are applied along the edges of intersected elements.

We will give details of the implementation of our method in FEniCS and will demonstrate the capability of our method on the rising bubble benchmark [4].

References [1] Burman, E., Claus, S., and Massing, A. (2015). A Stabilized Cut Finite Element Method for the Three Field Stokes Problem, SIAM J. Sci. Comput., Vol. 37(4), A1705-A1726. [2] Burman, E., Claus, S., Hansbo, P., Larson, M. G., and Massing, A. (2014). CutFEM: discretizing geometry and partial differential equations. International Journal for Numerical Methods in Engineering. 104(7), 472-501. [3] E. Burman (2010). Ghost penalty. Comptes Rendus Mathematique, Vol. 348 (21-22), 1217-1220. [4] Hysing, S. et al. (2009). Quantitative benchmark computations of two-dimensional bubble dynamics, International Journal for Numerical Methods in Fluids, Vol. 60 (11), 1259-1288.