120 research outputs found
HINT: Hierarchical Invertible Neural Transport for Density Estimation and Bayesian Inference
A large proportion of recent invertible neural architectures is based on a
coupling block design. It operates by dividing incoming variables into two
sub-spaces, one of which parameterizes an easily invertible (usually affine)
transformation that is applied to the other. While the Jacobian of such a
transformation is triangular, it is very sparse and thus may lack
expressiveness. This work presents a simple remedy by noting that (affine)
coupling can be repeated recursively within the resulting sub-spaces, leading
to an efficiently invertible block with dense triangular Jacobian. By
formulating our recursive coupling scheme via a hierarchical architecture, HINT
allows sampling from a joint distribution p(y,x) and the corresponding
posterior p(x|y) using a single invertible network. We demonstrate the power of
our method for density estimation and Bayesian inference on a novel data set of
2D shapes in Fourier parameterization, which enables consistent visualization
of samples for different dimensionalities
Matrix-free GPU implementation of a preconditioned conjugate gradient solver for anisotropic elliptic PDEs
Many problems in geophysical and atmospheric modelling require the fast
solution of elliptic partial differential equations (PDEs) in "flat" three
dimensional geometries. In particular, an anisotropic elliptic PDE for the
pressure correction has to be solved at every time step in the dynamical core
of many numerical weather prediction models, and equations of a very similar
structure arise in global ocean models, subsurface flow simulations and gas and
oil reservoir modelling. The elliptic solve is often the bottleneck of the
forecast, and an algorithmically optimal method has to be used and implemented
efficiently. Graphics Processing Units have been shown to be highly efficient
for a wide range of applications in scientific computing, and recently
iterative solvers have been parallelised on these architectures. We describe
the GPU implementation and optimisation of a Preconditioned Conjugate Gradient
(PCG) algorithm for the solution of a three dimensional anisotropic elliptic
PDE for the pressure correction in NWP. Our implementation exploits the strong
vertical anisotropy of the elliptic operator in the construction of a suitable
preconditioner. As the algorithm is memory bound, performance can be improved
significantly by reducing the amount of global memory access. We achieve this
by using a matrix-free implementation which does not require explicit storage
of the matrix and instead recalculates the local stencil. Global memory access
can also be reduced by rewriting the algorithm using loop fusion and we show
that this further reduces the runtime on the GPU. We demonstrate the
performance of our matrix-free GPU code by comparing it to a sequential CPU
implementation and to a matrix-explicit GPU code which uses existing libraries.
The absolute performance of the algorithm for different problem sizes is
quantified in terms of floating point throughput and global memory bandwidth.Comment: 18 pages, 7 figure
A hybrid Alternating Least Squares - TT Cross algorithm for parametric PDEs
We consider the approximate solution of parametric PDEs using the low-rank Tensor Train (TT) decomposition. Such parametric PDEs arise for example in uncertainty quantification problems in engineering applications. We propose an algorithm that is a hybrid of the alternating least squares and the TT cross methods. It computes a TT approximation of the whole solution, which is beneficial when multiple quantities of interest are sought. This might be needed, for example, for the computation of the probability density function (PDF) via the maximum entropy method [Kavehrad and Joseph, IEEE Trans. Comm., 1986]. The new algorithm exploits and preserves the block diagonal structure of the discretized operator in stochastic collocation schemes. This disentangles computations of the spatial and parametric degrees of freedom in the TT representation. In particular, it only requires solving independent PDEs at a few parameter values, thus allowing the use of existing high performance PDE solvers. In our numerical experiments, we apply the new algorithm to the stochastic diffusion equation and compare it with preconditioned steepest descent in the TT format, as well as with (multilevel) quasi-Monte Carlo and dimension-adaptive sparse grids methods. For sufficiently smooth random fields the new approach is orders of magnitude faster
Multilevel quasi-Monte Carlo for random elliptic eigenvalue problems I: Regularity and error analysis
Random eigenvalue problems are useful models for quantifying the uncertainty
in several applications from the physical sciences and engineering, e.g.,
structural vibration analysis, the criticality of a nuclear reactor or photonic
crystal structures. In this paper we present a simple multilevel quasi-Monte
Carlo (MLQMC) method for approximating the expectation of the minimal
eigenvalue of an elliptic eigenvalue problem with coefficients that are given
as a series expansion of countably-many stochastic parameters. The MLQMC
algorithm is based on a hierarchy of discretisations of the spatial domain and
truncations of the dimension of the stochastic parameter domain. To approximate
the expectations, randomly shifted lattice rules are employed. This paper is
primarily dedicated to giving a rigorous analysis of the error of this
algorithm. A key step in the error analysis requires bounds on the mixed
derivatives of the eigenfunction with respect to both the stochastic and
spatial variables simultaneously. An accompanying paper [Gilbert and Scheichl,
2021], focusses on practical extensions of the MLQMC algorithm to improve
efficiency, and presents numerical results
A Stein variational Newton method
Stein variational gradient descent (SVGD) was recently proposed as a general
purpose nonparametric variational inference algorithm [Liu & Wang, NIPS 2016]:
it minimizes the Kullback-Leibler divergence between the target distribution
and its approximation by implementing a form of functional gradient descent on
a reproducing kernel Hilbert space. In this paper, we accelerate and generalize
the SVGD algorithm by including second-order information, thereby approximating
a Newton-like iteration in function space. We also show how second-order
information can lead to more effective choices of kernel. We observe
significant computational gains over the original SVGD algorithm in multiple
test cases.Comment: 18 pages, 7 figure
- …