Programme of the 13rd Euro AD Workshop

Adjoint solvers are an essential ingredient in obtaining sensitivity information in Computational Fluid Dynamics. These gradients can then be used in numerical shape and topology optimisation, adaptive mesh refinement, flow control, uncertainty analysis and data assimilation,among many applications. AD is of course an eminently useful tool to derive adjoint codes for CFD, the research in this area has seen many important contributions. The presentation will report on the experience of the group in using Tapenade to derive adjoints for 3 particular codes. You will hear clear success stories which should encourage you to try this yourself on your limited-size in-house code where you can adapt the source code to perfect the collaboration with the AD tool. But you can also hear about our extended labour to upscale the same approach to a legacy CFD code where the size, complexity and legacy nature of the code makes this process much more difficult to achieve. The presentation will conclude with a discussion of a list of additional functionality that should be provided by AD tools in order to simplify the development and make-build process for legacy adjoint codes.

We discuss the application of Algorithmic or Automatic Differentiation (AD) to an industrial CFD code for shape optimization in a general context. AD is a method to generate discrete adjoint solvers in a semi-automated fashion. In an industrial framework, the differentiation process needs to be carefully designed to meet certain requirements on computer resources, coding standards and maintenance. Once a first adjoint solver is built in a semi-automated fashion, new developed code parts should be capable to be differentiated automatically and easy to be augmented to the adjoint solver.

TBA

Chebfun2 is a newly introduced extension to Chebfun, aimed at allowing users to compute numerically with functions of two variables. This talk describes first experiments of using AD for computing Fréchet derivatives of partial differential operators in two dimensions in Chebfun2, for solving nonlinear PDEs using Newton iteration and Chebyshev spectral methods. These experiments are inspired by the success such ideas have demonstrated for Chebfun's current capabilities of solving nonlinear ODEs.

GPU computing has become extremely popular in the recent years due to the sinking costs and the exploding efficiency of GPUs. The Vector-Forward mode of AD is highly suitable for computations on the GPU. The various direction vectors may be distributed among the cores and parallelizing the computation. GPU computation may also be used for evaluating the function and its directional derivatives at various points in the domain parallely. ADOL-C has been recently extended to be able to perform traceless forward computations on GPUs supporting the CUDA architecture. This talk will give details on the implementation and usage of the CUDA extension of ADOL-C for evaluating directional derivatives in parallel.

The MAD package (Forth, ACM Trans. Math. Softw. 32(2) 2006) facilitates overloaded forward mode AD in Matlab. Derivative combination operations are performed using high level array operations via the optimized derivvec class which ensures good performance for moderate problem size n. MAD was shown to outperform ADMAT (Verma, in Object Oriented Methods for Interoperable Scientific and Engineering Computing, SIAM, 1999). Kharche (PhD, Cranfield University, 2012) used source transformation to inline and specialise MAD's forward mode and derivative combination operations. For many problems his MSAD package demonstrated excellent performance. However, for the Burger's stiff ODE test problem and large n performance was some 30 times worse than compressed finite-differences. Recently Patterson et al (ACM Trans. Math. Softw, 39(3), 2013) have implemented forward ordered vertex-elimination AD using overloading to generate an evaluation trace which is then transformed using forward-ordered vertex elimination and used to write Matlab source code. On a small number of test cases they demonstrated performance in excess of that attained by MAD or ADMAT and, additionally, showed ADMAT outperforming MAD for large n. In this presentation we will briefly review MAD's forward mode and the approach of MSAD. We will then show some performance statistics that suggest some routes to performance improvements for the derivvec class. In order to ensure software quality we outline our new unit testing suite that makes use of the unit testing facilities of Matlab Release 2013a. We hope to present improved performance statistics!

We discuss the structure of nonlinear functions and mechanisms to exploit this structure. We describe progress in the automatic detection and exploitation of partial separability. We also examine the nonlinearity structure in power grid analysis and discretized partial differential equations.

The FEniCS project is a collection of free software for the automated solution of partial differential equations. Users specify the variational formulation of their finite element discretisation in an expressive domain-specific language, which is used by a finite element form compiler to automatically generate the code for its assembly. Uniquely, this representation retains the mathematical structure of a finite element model in a symbolically manipulable form, greatly facilitating the derivation of its associated tangent linear and adjoint models. This naturally fits with a high-level approach to algorithmic differentiation, where individual linear and nonlinear solves are included in the set of elementary operations. In this talk we discuss the high-level approach to algorithmic differentiation taken, and demonstrate its advantages. In particular, the approach is almost fully automatic, achieves excellent performance, and works naturally in parallel.

We show that all piecewise-linear functions in n variables can be brought into abs-normal form with a maximal switching depth of 2n-1. On the abs-normal form we can perform some standard numerical analysis based on classical linear algebra. Also PL systems of equations turn into linear complementarity problems.

The prevailing wisdom is that it is not feasible to differentiate naively though equation solution algorithms that have a singularity at their solution. We present an analysis which suggests that this may be possible after all. The key issue turns out to be deciding when a coefficient is "actually" zero: is the real part of (a + bt)/(c + dt) equal to a/c, or to b/d?

TBA

Computer aided simulation of guided acoustic waves in single- or multilayered waveguides is an essential tool for several applications of acoustics and ultrasonics (i.e. pipe inspection, noise reduction). To simulate wave propagation in geometrically simple waveguides (plates or rods), one may employ the analytical global matrix method. This requires the computation of all roots of the determinate of a certain submatrix. The evaluation of all real or even complex roots is actually the methods most concerning restriction. Previous approaches base on so called mode-tracers which use the physical phenomenon that solutions (roots) appear in a certain pattern (waveguide modes) and thus use known solutions to limit the root finding algorithms searchspace with respect to consecutive solutions. As the limitation of searchspace might be unstable in some cases, we propose to replace the mode-tracer with a suitable version of an interval Newton method based on Intlab. To apply this interval based method, we extended the interval and derivative computation provided by Intlab such that corresponding information is also available for Bessel functions used in the circular model (rods) of acoustic waveguides. We present numerical results of a simple acoustic waveguide and discuss extensions required for more realistic scenarios.

Least squares regression is among the most frequently used method for estimating the parameters of financial models. It is necessary to find the global optima of these functions in order to successfully meet the least squares. Unfortunately, in case of ill-conditioned functions, most existing methods often converge to local minimum. This presentation investigates a hybrid approach for non-linear least squares global optimization based on Interval Analysis and Automatic Differentiation methods. Additional heuristic algorithms and local optimization solvers are considered to narrow the input domain. The proposed algorithm is applied in the calibration of Heston model for option pricing (a semi-closed form) to fit market quotes. Furthermore, the parallel approach using GPU is briefly discussed.

The computations of the high-order partial derivatives in a given problem are in general tedious or not accurate. To combat such shortcomings, a new method for calculating exact high-order sensitivities using multi-complex numbers is presented. Inspired by the recent complex step method that is only valid for first order sensitivities, the new multi-complex approach is valid to arbitrary order. The mathematical theory behind this approach is revealed, and an efficient procedure for the automatic implementation of the method is described. Several applications are presented to validate and demonstrate the accuracy and efficiency of the algorithm. The results are compared to conventional approaches such as finite differencing, the complex step method, and two separate automatic differentiation tools. Our multi-complex method is shown to have many advantages, and it is therefore expected to be useful for any algorithm exploiting high-order derivatives, such as many non-linear programming solvers.

MPI is the de facto standard in distributed parallel computing. In particular, its implementations are used to parallize numerical simulations on computer clusters. Any application of AD on such a simulation code will probably have to deal with the MPI communication.

In this talk, a short list Dos and Don'ts is presented with regard to adjoining MPI. Furthermore, we discuss the current state of a generic adjoint MPI library.

The adjoint mode implies a reversal of the data dependencies and consequently a reversal of communications in parallelized models. Building on previous studies regarding the adjoining of MPI two-sided communications, we investigate the construction of adjoints for certain one-sided MPI communications and the overhead their use implies for the adjoint.