ABSTRACT. In the first part of the talk I will describe in general terms the link between classical optimisation techniques and the Bayesian approach to statistical inversion as outlined in the seminal book of [Kaipio and Somersalo, 2005]. Under the assumption of an additive Gaussian noise model, a Gaussian prior distribution and a linear parameter-to-observable map, it is possible to uniquely characterise the Bayesian posterior as Gaussian with the maximum aposteriori (MAP) point equal to the minimum of a classic regularised minimisation problem and covariance matrix equal to the inverse of the Hessian of the functional evaluated at the MAP point. I will also discuss techniques that can be used when these assumptions break down.

In the second part of the talk I will describe a method implemented within dolfin-adjoint [Funke and Farrell, arXiv 2013] to quantify the uncertainty in the recovered material parameters of a hyperelastic solid from partial and noisy observations of the displacement field in the domain. The finite element discretisation of the adjoint and higher-order adjoint (Hessian) equations are derived automatically from the high-level UFL representation of the problem. The resulting equations are solved using PETSc. I will concentrate on finding the eigenvalue decomposition of the posterior covariance matrix (Hessian). The eigenvectors associated with the lowest eigenvalues of the Hessian correspond with the directions in parameter space least constrained by the observations [Flath et al. 2011]. This eigenvalue problem is tricky to solve efficiently because the Hessian is very large (on the order of the number of parameters) and dense (meaning that only its action on a vector can be calculated, each involving considerable expense).

Finally, I will show some illustrative examples including the uncertainty associated with deriving the material properties of a 3D hyperelastic block with a stiff inclusion with knowledge only of the displacements on the boundary of the domain.

J. Kaipio and E. Somersalo, Statistical and Computational Inverse Problems, vol. 160. New York: Springer-Verlag, 2005. S. W. Funke and P. E. Farrell, “A framework for automated PDE-constrained optimisation,” arXiv:1302.3894 [cs], Feb. 2013. H. P. Flath, L. C. Wilcox, V. Akçelik, J. Hill, B. van Bloemen Waanders, and O. Ghattas, “Fast Algorithms for Bayesian Uncertainty Quantification in Large-Scale Linear Inverse Problems Based on Low-Rank Partial Hessian Approximations,” SIAM J. Sci. Comput., vol. 33, no. 1, pp. 407–432, Feb. 2011.

ABSTRACT. Non-destructive testing and evaluation (NDT&E) by ultrasonic (or electromagnetic) waves is very important for the identification of material parameters, for quality assurance and life-time control. New materials or composites, complicated geometries and higher resolution requirements trigger a need for improved NDT&E techniques and procedures. A priori knowledge and expertise from prototyping  are generally expensive or even unavailable. Numerical field simulation combined with adjoint-based optimizing techniques are a flexible alternative.

In this presentation, I will discuss the challenges involved with this problem and my progress in using and adapting FEniCS and dolfin-adjoint for overcoming these challenges. The main issues that were overcome so far are solving Discontinuous Galerkin systems more efficiently and time-dependent objective functionals that use reference data in certain points in the considered geometry.

ABSTRACT. We consider a saddle point formulation of a system consisting of an elliptic PDE on a two- or three-dimensional domain coupled to a differential equation on a domain of lower topological dimension. To solve the problem with FEniCS, the missing support for mixed function spaces with components defined over meshes of different topological dimension must be addressed. A memory efficient alternative to constraining the excess degrees of freedom is the explicitly constructed trace operator to the lower-dimensional domain. The constructed trace operator enables us to use FEniCS to explore various spaces/inner products in which the constraint is enforced. It is found that with the natural $L^2$ inner product the efficient preconditioners for the resulting linear systems are based on fractional Sobolev norms.

ABSTRACT. In many applications the geometry of the domain is characterized by one dimension being much smaller than the other two. After appropriate expansions in orders of an aspect ratio are truncated one is usually left with a two dimensional model involving the large dimensions only. In such cases the manifold technology within FEniCS becomes very useful. In many real world applications however the resulting two dimensional domain consists of several, and in some cases many hundreds or thousands, of intersecting surfaces. Such domains are not manifolds, but we will here demonstrate how the manifolds machinery in FEniCS allows us to define and solve variational problems on these domains. We demonstrate the solution of an elliptic equation as well as a time dependent, nonlinear advection-diffusion equation on several examples of such a domain.

ABSTRACT. Shell, plate and beam (thin) structures are widely used in civil, mechanical and aeronautical engineering because they are capable of carrying high loads with a minimal amount of structural mass. Because the out-of-plane dimension is usually much smaller than the two in-plane dimensions, it is possible to asymptotically reduce the full 3D equations of elasticity to a whole variety of equivalent 2D models posed on a manifold embedded in 3D space. This reduction results in massively reduced computational expense and remains a necessity for practical large-scale computation of structures of real engineering interest such as the fuselage of an aircraft.

The numerical solution of such mathematical models is a challenging task, especially for very thin shells when shear and membrane locking effects require special attention. As originally noted by [Hale and Baiz, 2013], the high-level form language UFL provides an excellent framework for writing extensible, reusable and pedagogical numerical models of thin structures. To our knowledge fenics-shells represents the first unified open-source implementation of a wide range of thin structural models, including Reissner-Mindlin, Kirchhoff-Love, Von Karman and hierarchical (higher-order) plates, and Madare-Naghdi and Madare-Koiter shell models.

Because of the broad scope of fenics-shells, in this talk we will focus on how to cure numerical locking by applying the Mixed Interpolation of Tensorial Components (MITC) approach of [Dvorkin and Bathe, 1986] and [Lee and Bathe, 2010] to a shell with an initially flat reference configuration. The MITC approach consists of an element-by-element interpolation of the degrees of freedom of the rotations onto the degrees of freedom of a reduced rotation space, the latter typically constructed using H(curl) conforming finite elements such as the rotated Raviart-Thomas-Nédélec elements. Then, the bilinear form is constructed on the underlying H(curl) space. Because of the interpolation operator, the original problem is expressed in terms of the degrees of freedom for the rotations only. Within DOLFIN we have implemented this projection operation using two UFL forms within a custom assembler compiled just-in-time using Instant. We show numerical convergence studies that match the apriori bounds available in the literature.

E. N. Dvorkin and K.-J. Bathe, “A continuum mechanics based four-node shell element for general non-linear analysis,” Engineering Computations, vol. 1, no. 1, pp. 77–88, 1984. P. S. Lee and K. J. Bathe, “The quadratic MITC plate and MITC shell elements in plate bending,” Advances in Engineering Software, vol. 41, no. 5, pp. 712–728, 2010. J. S. Hale and P. M. Baiz, “Towards effective shell modelling with the FEniCS project” presented at the FEniCS Conference 2013, Jesus College, Cambridge, 19-Mar-2013.

ABSTRACT. The Large Eddy Simulation (LES) of turbulent flows in complex geometries remains a key challenge for engineering applications. While LES has historically relied on physical modelling of the small scales, the approach denoted as Implicit LES, consists in seeing the turbulence model as a by-product of the discretization of the Navier-Stokes or Euler equations.

Reference methods [1] usually rely on an equal order Lagrange finite element discretization of the Navier–Stokes equation stabilized by a Streamline-Diffusion/SUPG type term. Nonlinear viscosities based on the entropy residual [2], motivated by the theory of measure-valued solution for conservation laws, are also considered. In this case, Navier–Stokes equations are seen as a elliptic regularization of the Euler equations: being based on the sole computation of a residual-based ``eddy viscosity'', they are convenient to insert in splitting and explicit schemes [3].

Applications to incompressible turbulent flows, in the homogeneous and variable density case, as well as compressible flows are presented together with the recent developments in DOLFIN-HPC [5] that they motivated.

[1] J. Hoffman, C. Johnson. A new approach to Computational Turbulence Modeling. Comput. Methods Appl. Mech. Engrg., 195, pp2865–2880, (2006).

[2] J. L. Guermond, R. Pasquetti, B. Popov. From suitable weak solutions to entropy viscosity. J. Scientific Comput., 49, pp35–50, (2011).

[3] M. Nazarov, J. Hoffman. Residual based artificial viscosity for simulation of turbulent compressible flow using adaptive finite element methods. Int. J. Num. Methods Fluids, 71(3): 339-357, (2013).

[4] A. Larcher, K. Müller, M. Nazarov, J. Hoffman. A residual viscosity-based splitting scheme for variable density flows, KTH preprint, (2015).

[5] N. Jansson, J. Hoffman, J. Jansson. Framework for Massively Par- allel Adaptive Finite Element CFD on Tetrahedral Meshes. SIAM J. Sci. Comput., 34, (2012).

ABSTRACT. We present work where FEniCS is used to simulate the blood flow in the left ventricle of the human heart. A workflow is established where the heart model is calibrated to patient-specific data from ultrasound measurements, to enable individualised simulations. An Arbitrary Lagrangian-Eulerian (ALE) finite element method is used to model the blood flow, which is implemented as a Unicorn solver based on Dolfin-HPC. A fluid-structure interaction (FSI) model of the aortic valves is also developed, which can be used to model both biological and mechanical valves. The goal is to develop a clinical decision support tool to simulate heart disease progression, and predictive outcomes of treatment.

ABSTRACT. When modelling carbon dioxide (CO2) sequestration, the process of storing CO2 underground in layers of porous rock as a way of mitigating its negative impact on the environment, it is necessary to solve a nonlinear system of partial differential equations. FEniCS has been used to simulate a model system of CO2 storage in the presence of background flow in an aquifer. For this system, in which the background flow is in the opposite direction to the diffusive flux, the choice of basis functions is critical for creating a stable solution method. Moreover, there is a compatibility condition between the function space for the CO2 concentration and the velocity space to avoid spurious density driven velocity fluctuations. FEniCS makes changing the basis functions for different fields straightforward and rapid, which in turn permits different solution strategies to be explored. Solutions computed with FEniCS are compared to a steady, one-dimensional analytical model, and validate the analytical model for long-term CO2 transport against a background flow.

ABSTRACT. A novel class of preconditioners in the PETSc library that are based on Balancing Domain Decomposition by Constraints (BDDC) methods is presented. The talk introduces the BDDC algorithm and its current implementation in PETSc, providing details of available user customizations. An experimental interface to the Finite Element Tearing and Interconnecting Dual-Primal (FETI-DP) method in PETSc is also discussed. Numerical results for Raviart-Thomas and Nedelec finite elements are provided, showing robustness of the methods with respect to heterogeneous distributions of the coefficients of elliptic PDEs.

ABSTRACT. For a wide class of problems, both linear and nonlinear, multigrid is an optimal solution method. Along with black and grey-box algebraic approaches, when the underlying mesh has some (refinement) structure, geometric multigrid can be highly effective. As well as the prospect of solving linear systems faster (and with a smaller memory footprint) than AMG, geometric multigrid can be applied directly to nonlinear problems, negating the need for global linearisation (à la Newton).

In this talk, I will present recent progress in providing a framework for developing GMG solvers in Firedrake. As well as providing the building blocks to write such solvers "by hand", we exploit PETSc's multigrid framework to provide "plug and play" multigrid if a mesh hierarchy is available. I will describe how the composition of PETSc's linear algebra abstractions and the symbolic information available in UFL make this possible, illustrating with some examples.