BZPEERPapers, people, evidence and research paths
Menu

Discover

SearchFind papers, people, topics, institutions and opportunities from one place.GraphSee how papers, authors, topics and evidence connect.

Research tools

HomeContinue from the papers, people and topics that matter now.CollectionsTurn saved papers into reading lists, evidence maps and shareable context.ResearchReconstruct the path that led you from search to paper to question.CollaborationMove from a paper, topic or expert request into a focused thread.PublishCreate full-text research articles with metadata, formulas and versions.

Network

ExpertsFind specialists whose work and reviews explain why they can help.CommunitiesJoin research groups organized around papers, topics and shared evidence.

Opportunities

ListingsBrowse roles, calls and collaborations connected to real research context.

Account

Log inReturn to your saved papers, graph views and active threads.Create accountStart saving papers, building a research trail and joining teams.
Log inCreate account

Source author record

Robert Scheichl

Robert Scheichl appears in the imported research catalog. Authorship, coauthor and topic links are available while profile ownership is still unclaimed.

ResearcherUnclaimed source record
math.NANumerical AnalysisDistributed, Parallel, and Cluster Computingmath.PRMathematical SoftwareComputationComputational Engineering, Finance, and Sciencecond-mat.softcond-mat.stat-mechhep-latmath.STMethodologyphysics.ao-phphysics.comp-phStatistics Theory

Catalog footprint

What is connected

19works
15topics
4close collaborators

Actions

Connect this record

Log in to claim
Open graphBrowse works

Research graph

See the researcher in context

Open full explorer

Inspect adjacent papers, topics, institutions and collaborators without losing the researcher page.

Building this map preview

BZPEER is loading the nearby papers, people, topics and institutions for this page.

Published work

19 published item(s)

preprint2022arXiv

Multilevel Delayed Acceptance MCMC

We develop a novel Markov chain Monte Carlo (MCMC) method that exploits a hierarchy of models of increasing complexity to efficiently generate samples from an unnormalized target distribution. Broadly, the method rewrites the Multilevel MCMC approach of Dodwell et al. (2015) in terms of the Delayed Acceptance (DA) MCMC of Christen & Fox (2005). In particular, DA is extended to use a hierarchy of models of arbitrary depth, and allow subchains of arbitrary length. We show that the algorithm satisfies detailed balance, hence is ergodic for the target distribution. Furthermore, multilevel variance reduction is derived that exploits the multiple levels and subchains, and an adaptive multilevel correction to coarse-level biases is developed. Three numerical examples of Bayesian inverse problems are presented that demonstrate the advantages of these novel methods. The software and examples are available in PyMC3.

0 citations0 reviewsNo code linkedOpen access
preprint2022arXiv

Wavenumber explicit convergence of a multiscale generalized finite element method for heterogeneous Helmholtz problems

In this paper, a generalized finite element method (GFEM) with optimal local approximation spaces for solving high-frequency heterogeneous Helmholtz problems is systematically studied. The local spaces are built from selected eigenvectors of carefully designed local eigenvalue problems defined on generalized harmonic spaces. At both continuous and discrete levels, $(i)$ wavenumber explicit and nearly exponential decay rates for local and global approximation errors are obtained without any assumption on the size of subdomains; $(ii)$ a quasi-optimal convergence of the method is established by assuming that the size of subdomains is $O(1/k)$ ($k$ is the wavenumber). A novel resonance effect between the wavenumber and the dimension of local spaces on the decay of error with respect to the oversampling size is implied by the analysis. Furthermore, for fixed dimensions of local spaces, the discrete local errors are proved to converge as $h\rightarrow 0$ ($h$ denoting the mesh size) towards the continuous local errors. The method at the continuous level extends the plane wave partition of unity method [I. Babuska and J. M. Melenk, Int.\;J.\;Numer.\;Methods Eng., 40 (1997), pp.~727--758] to the heterogeneous-coefficients case, and at the discrete level, it delivers an efficient non-iterative domain decomposition method for solving discrete Helmholtz problems resulting from standard FE discretizations. Numerical results are provided to confirm the theoretical analysis and to validate the proposed method.

0 citations0 reviewsNo code linkedOpen access
preprint2020arXiv

A High-Performance Implementation of a Robust Preconditioner for Heterogeneous Problems

We present an efficient implementation of the highly robust and scalable GenEO preconditioner in the high-performance PDE framework DUNE. The GenEO coarse space is constructed by combining low energy solutions of a local generalised eigenproblem using a partition of unity. In this paper we demonstrate both weak and strong scaling for the GenEO solver on over 15,000 cores by solving an industrially motivated problem with over 200 million degrees of freedom. Further, we show that for highly complex parameter distributions arising in certain real-world applications, established methods become intractable while GenEO remains fully effective. The purpose of this paper is two-fold: to demonstrate the robustness and high parallel efficiency of the solver and to document the technical details that are crucial to the efficiency of the code.

0 citations0 reviewsNo code linkedOpen access
preprint2020arXiv

Full Error Analysis and Uncertainty Quantification for the Heterogeneous Transport Equation in Slab Geometry

We present an analysis of multilevel Monte Carlo techniques for the forward problem of uncertainty quantification for the radiative transport equation, when the coefficients ({\em cross-sections}) are heterogenous random fields. To do this, we first give a new error analysis for the combined spatial and angular discretisation in the deterministic case, with error estimates which are explicit in the coefficients (and allow for very low regularity and jumps). This detailed error analysis is done for the 1D space - 1D angle slab geometry case with classical diamond differencing. Under reasonable assumptions on the statistics of the coefficients, we then prove an error estimate for the random problem in a suitable Bochner space. Because the problem is not self-adjoint, stability can only be proved under a path-dependent mesh resolution condition. This means that, while the Bochner space error estimate is of order $\mathcal{O}(h^η)$ for some $η$, where $h$ is a (deterministically chosen) mesh diameter, smaller mesh sizes might be needed for some realisations. Under reasonable assumptions we show that the expected cost for computing a typical quantity of interest remains of the same order as for a single sample. This leads to rigorous complexity estimates for Monte Carlo and multilevel Monte Carlo: For particular linear solvers, the multilevel version gives up to two orders of magnitude improvement over Monte Carlo. We provide numerical results supporting the theory.

0 citations0 reviewsNo code linkedOpen access
preprint2020arXiv

Multilevel Monte Carlo for quantum mechanics on a lattice

Monte Carlo simulations of quantum field theories on a lattice become increasingly expensive as the continuum limit is approached since the cost per independent sample grows with a high power of the inverse lattice spacing. Simulations on fine lattices suffer from critical slowdown, the rapid growth of autocorrelations in the Markov chain. This causes a strong increase in the number of lattice configurations that have to be generated to obtain statistically significant results. This paper discusses hierarchical sampling methods to tame the growth in autocorrelations. Combined with multilevel variance reduction, this significantly reduces the computational cost of simulations for given tolerances $ε_{\text{disc}}$ on the discretisation error and $ε_{\text{stat}}$ on the statistical error. For observables with lattice errors of order $α$ and integrated autocorrelation times that grow like $τ_{\mathrm{int}}\propto a^{-z}$, multilevel Monte Carlo (MLMC) reduces the cost from $\mathcal{O}(ε_{\text{stat}}^{-2}ε_{\text{disc}}^{-(1+z)/α})$ to $\mathcal{O}(ε_{\text{stat}}^{-2}\vert\log ε_{\text{disc}} \vert^2+ε_{\text{disc}}^{-1/α})$ or $\mathcal{O}(ε_{\text{stat}}^{-2}+ε_{\text{disc}}^{-1/α})$. Higher gains are expected for simulations of quantum field theories in $D$ dimensions. The efficiency of the approach is demonstrated on two model systems, including a topological oscillator that is badly affected by critical slowdown from topological charge freezing. On fine lattices, the new methods are orders of magnitude faster than standard Hybrid Monte Carlo sampling. For high resolutions, MLMC can be used to accelerate even the cluster algorithm for the topological oscillator. Performance is further improved through perturbative matching which guarantees efficient coupling of theories on the multilevel hierarchy.

0 citations0 reviewsNo code linkedOpen access
preprint2020arXiv

Unified Analysis of Periodization-Based Sampling Methods for Matérn Covariances

The periodization of a stationary Gaussian random field on a sufficiently large torus comprising the spatial domain of interest is the basis of various efficient computational methods, such as the classical circulant embedding technique using the fast Fourier transform for generating samples on uniform grids. For the family of Matérn covariances with smoothness index $ν$ and correlation length $λ$, we analyse the nonsmooth periodization (corresponding to classical circulant embedding) and an alternative procedure using a smooth truncation of the covariance function. We solve two open problems: the first concerning the $ν$-dependent asymptotic decay of eigenvalues of the resulting circulant in the nonsmooth case, the second concerning the required size in terms of $ν$, $λ$ of the torus when using a smooth periodization. In doing this we arrive at a complete characterisation of the performance of these two approaches. Both our theoretical estimates and the numerical tests provided here show substantial advantages of smooth truncation.

0 citations0 reviewsNo code linkedOpen access
preprint2019arXiv

A fully adaptive multilevel stochastic collocation strategy for solving elliptic PDEs with random data

We propose and analyse a fully adaptive strategy for solving elliptic PDEs with random data in this work. A hierarchical sequence of adaptive mesh refinements for the spatial approximation is combined with adaptive anisotropic sparse Smolyak grids in the stochastic space in such a way as to minimize the computational cost. The novel aspect of our strategy is that the hierarchy of spatial approximations is sample dependent so that the computational effort at each collocation point can be optimised individually. We outline a rigorous analysis for the convergence and computational complexity of the adaptive multilevel algorithm and we provide optimal choices for error tolerances at each level. Two numerical examples demonstrate the reliability of the error control and the significant decrease in the complexity that arises when compared to single level algorithms and multilevel algorithms that employ adaptivity solely in the spatial discretisation or in the collocation procedure.

0 citations0 reviewsNo code linkedOpen access
preprint2019arXiv

Correction of coarse-graining errors by a two-level method: application to the Asakura-Oosawa model

We present a method that exploits self-consistent simulation of coarse-grained and fine-grained models, in order to analyse properties of physical systems. The method uses the coarse-grained model to obtain a first estimate of the quantity of interest, before computing a correction by analysing properties of the fine system. We illustrate the method by applying it to the Asakura-Oosawa (AO) model of colloid-polymer mixtures. We show that the liquid-vapour critical point in that system is affected by three-body interactions which are neglected in the corresponding coarse-grained model. We analyse the size of this effect and the nature of the three-body interactions. We also analyse the accuracy of the method, as a function of the associated computational effort.

0 citations0 reviewsNo code linkedOpen access
preprint2019arXiv

High-performance dune modules for solving large-scale, strongly anisotropic elliptic problems with applications to aerospace composites

The key innovation in this paper is an open-source, high-performance iterative solver for high contrast, strongly anisotropic elliptic partial differential equations implemented within dune-pdelab. The iterative solver exploits a robust, scalable two-level additive Schwarz preconditioner, GenEO (Spillane et al. 2014). The development of this solver has been motivated by the need to overcome the limitations of commercially available modeling tools for solving structural analysis simulations in aerospace composite applications. Our software toolbox dune-composites encapsulates the mathematical complexities of the underlying packages within an efficient C++ framework, providing an application interface to our new high-performance solver. We illustrate its use on a range of industrially motivated examples, which should enable other scientists to build on and extend dune-composites and the GenEO preconditioner for use in their own applications. We demonstrate the scalability of the solver on more than 15,000 cores of the UK national supercomputer Archer, solving an aerospace composite problem with over 200 million degrees of freedom in a few minutes. This scale of computation brings composites problems that would otherwise be unthinkable into the feasible range. To demonstrate the wider applicability of the new solver, we also confirm the robustness and scalability of the solver on SPE10, a challenging benchmark in subsurface flow/reservoir simulation.

0 citations0 reviewsNo code linkedOpen access
preprint2016arXiv

Multilevel Quasi-Monte Carlo Methods for Lognormal Diffusion Problems

In this paper we present a rigorous cost and error analysis of a multilevel estimator based on randomly shifted Quasi-Monte Carlo (QMC) lattice rules for lognormal diffusion problems. These problems are motivated by uncertainty quantification problems in subsurface flow. We extend the convergence analysis in [Graham et al., Numer. Math. 2014] to multilevel Quasi-Monte Carlo finite element discretizations and give a constructive proof of the dimension-independent convergence of the QMC rules. More precisely, we provide suitable parameters for the construction of such rules that yield the required variance reduction for the multilevel scheme to achieve an $\varepsilon$-error with a cost of $\mathcal{O}(\varepsilon^{-θ})$ with $θ< 2$, and in practice even $θ\approx 1$, for sufficiently fast decaying covariance kernels of the underlying Gaussian random field inputs. This confirms that the computational gains due to the application of multilevel sampling methods and the gains due to the application of QMC methods, both demonstrated in earlier works for the same model problem, are complementary. A series of numerical experiments confirms these gains. The results show that in practice the multilevel QMC method consistently outperforms both the multilevel MC method and the single-level variants even for non-smooth problems.

0 citations0 reviewsNo code linkedOpen access
preprint2016arXiv

Robust Numerical Upscaling of Elliptic Multiscale Problems at High Contrast

We present a new approach to the numerical upscaling for elliptic problems with rough diffusion coefficient at high contrast. It is based on the localizable orthogonal decomposition of $H^1$ into the image and the kernel of some novel stable quasi-interpolation operators with local $L^2$-approximation properties, independent of the contrast. We identify a set of sufficient assumptions on these quasi-interpolation operators that guarantee in principle optimal convergence without pre-asymptotic effects for high-contrast coefficients. We then give an example of a suitable operator and establish the assumptions for a particular class of high-contrast coefficients. So far this is not possible without any pre-asymptotic effects, but the optimal convergence is independent of the contrast and the asymptotic range is largely improved over other discretisation schemes. The new framework is sufficiently flexible to allow also for other choices of quasi-interpolation operators and the potential for fully robust numerical upscaling at high contrast.

0 citations0 reviewsNo code linkedOpen access
preprint2016arXiv

Scheduling massively parallel multigrid for multilevel Monte Carlo methods

The computational complexity of naive, sampling-based uncertainty quantification for 3D partial differential equations is extremely high. Multilevel approaches, such as multilevel Monte Carlo (MLMC), can reduce the complexity significantly, but to exploit them fully in a parallel environment, sophisticated scheduling strategies are needed. Often fast algorithms that are executed in parallel are essential to compute fine level samples in 3D, whereas to compute individual coarse level samples only moderate numbers of processors can be employed efficiently. We make use of multiple instances of a parallel multigrid solver combined with advanced load balancing techniques. In particular, we optimize the concurrent execution across the three layers of the MLMC method: parallelization across levels, across samples, and across the spatial grid. The overall efficiency and performance of these methods will be analyzed. Here the scalability window of the multigrid solver is revealed as being essential, i.e., the property that the solution can be computed with a range of process numbers while maintaining good parallel efficiency. We evaluate the new scheduling strategies in a series of numerical tests, and conclude the paper demonstrating large 3D scaling experiments.

0 citations0 reviewsNo code linkedOpen access
preprint2015arXiv

Efficient Multigrid Preconditioners for Atmospheric Flow Simulations at High Aspect Ratio

Many problems in fluid modelling require the efficient solution of highly anisotropic elliptic partial differential equations (PDEs) in "flat" domains. For example, in numerical weather- and climate-prediction an elliptic PDE for the pressure correction has to be solved at every time step in a thin spherical shell representing the global atmosphere. This elliptic solve can be one of the computationally most demanding components in semi-implicit semi-Lagrangian time stepping methods which are very popular as they allow for larger model time steps and better overall performance. With increasing model resolution, algorithmically efficient and scalable algorithms are essential to run the code under tight operational time constraints. We discuss the theory and practical application of bespoke geometric multigrid preconditioners for equations of this type. The algorithms deal with the strong anisotropy in the vertical direction by using the tensor-product approach originally analysed by Börm and Hiptmair [Numer. Algorithms, 26/3 (2001), pp. 219-234]. We extend the analysis to three dimensions under slightly weakened assumptions, and numerically demonstrate its efficiency for the solution of the elliptic PDE for the global pressure correction in atmospheric forecast models. For this we compare the performance of different multigrid preconditioners on a tensor-product grid with a semi-structured and quasi-uniform horizontal mesh and a one dimensional vertical grid. The code is implemented in the Distributed and Unified Numerics Environment (DUNE), which provides an easy-to-use and scalable environment for algorithms operating on tensor-product grids. Parallel scalability of our solvers on up to 20,480 cores is demonstrated on the HECToR supercomputer.

0 citations0 reviewsNo code linkedOpen access
preprint2015arXiv

Petascale elliptic solvers for anisotropic PDEs on GPU clusters

Memory bound applications such as solvers for large sparse systems of equations remain a challenge for GPUs. Fast solvers should be based on numerically efficient algorithms and implemented such that global memory access is minimised. To solve systems with up to one trillion ($10^{12}$) unknowns the code has to make efficient use of several million individual processor cores on large GPU clusters. We describe the multi-GPU implementation of two algorithmically optimal iterative solvers for anisotropic elliptic PDEs which are encountered in atmospheric modelling. In this application the condition number is large but independent of the grid resolution and both methods are asymptotically optimal, albeit with different absolute performance. We parallelise the solvers and adapt them to the specific features of GPU architectures, paying particular attention to efficient global memory access. We achieve a performance of up to 0.78 PFLOPs when solving an equation with $0.55\cdot 10^{12}$ unknowns on 16384 GPUs; this corresponds to about $3\%$ of the theoretical peak performance of the machine and we use more than $40\%$ of the peak memory bandwidth with a Conjugate Gradient (CG) solver. Although the other solver, a geometric multigrid algorithm, has a slightly worse performance in terms of FLOPs per second, overall it is faster as it needs less iterations to converge; the multigrid algorithm can solve a linear PDE with half a trillion unknowns in about one second.

0 citations0 reviewsNo code linkedOpen access
preprint2013arXiv

An Euler-Poisson Scheme for Lévy driven SDEs

We describe an Euler scheme to approximate solutions of Lévy driven Stochastic Differential Equations (SDE) where the grid points are random and given by the arrival times of a Poisson process. This result extends a previous work of the authors in Ferreiro-Castilla et al. (2012). We provide a complete numerical analysis of the algorithm to approximate the terminal value of the SDE and proof that the approximation converges in mean square error with rate $\mathcal{O}(n^{-1/2})$. The only requirement of the methodology is to have exact samples from the resolvent of the Lévy process driving the SDE; classic examples such as stable processes, subclasses of spectrally one sided Lévy processes and new families such as meromorphic Lévy processes (cf. Kuznetsov et al. (2011)) are some examples for which the implementation of our algorithm is straightforward.

0 citations0 reviewsNo code linkedOpen access
preprint2013arXiv

Massively parallel solvers for elliptic PDEs in Numerical Weather- and Climate Prediction

The demand for substantial increases in the spatial resolution of global weather- and climate- prediction models makes it necessary to use numerically efficient and highly scalable algorithms to solve the equations of large scale atmospheric fluid dynamics. For stability and efficiency reasons several of the operational forecasting centres, in particular the Met Office and the ECMWF in the UK, use semi-implicit semi-Lagrangian time stepping in the dynamical core of the model. The additional burden with this approach is that a three dimensional elliptic partial differential equation (PDE) for the pressure correction has to be solved at every model time step and this often constitutes a significant proportion of the time spent in the dynamical core. To run within tight operational time scales the solver has to be parallelised and there seems to be a (perceived) misconception that elliptic solvers do not scale to large processor counts and hence implicit time stepping can not be used in very high resolution global models. After reviewing several methods for solving the elliptic PDE for the pressure correction and their application in atmospheric models we demonstrate the performance and very good scalability of Krylov subspace solvers and multigrid algorithms for a representative model equation with more than $10^{10}$ unknowns on 65536 cores on HECToR, the UK's national supercomputer. For this we tested and optimised solvers from two existing numerical libraries (DUNE and hypre) and implemented both a Conjugate Gradient solver and a geometric multigrid algorithm based on a tensor-product approach which exploits the strong vertical anisotropy of the discretised equation. We study both weak and strong scalability and compare the absolute solution times for all methods; in contrast to one-level methods the multigrid solver is robust with respect to parameter variations.

0 citations0 reviewsNo code linkedOpen access
preprint2013arXiv

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.

0 citations0 reviewsNo code linkedOpen access
preprint2013arXiv

Mixed Finite Element Analysis of Lognormal Diffusion and Multilevel Monte Carlo Methods

This work is motivated by the need to develop efficient tools for uncertainty quantification in subsurface flows associated with radioactive waste disposal studies. We consider single phase flow problems in random porous media described by correlated lognormal distributions. We are interested in the error introduced by a finite element discretisation of these problems. In contrast to several recent works on the analysis of standard nodal finite element discretisations, we consider here mass-conservative lowest order Raviart-Thomas mixed finite elements. This is very important since local mass conservation is highly desirable in realistic groundwater flow problems. Due to the limited spatial regularity and the lack of uniform ellipticity and boundedness of the operator the analysis is non-trivial in the presence of lognormal random fields. We establish finite element error bounds for Darcy velocity and pressure, as well as for a more accurate recovered pressure approximation. We then apply the error bounds to prove convergence of the multilevel Monte Carlo algorithm for estimating statistics of these quantities. Moreover, we prove convergence for a class of bounded, linear functionals of the Darcy velocity. An important special case is the approximation of the effective permeability in a 2D flow cell. We perform numerical experiments to confirm the convergence results.

0 citations0 reviewsNo code linkedOpen access
preprint2013arXiv

Multilevel Monte Carlo simulation for Levy processes based on the Wiener-Hopf factorisation

In Kuznetsov et al. (2011) a new Monte Carlo simulation technique was introduced for a large family of Levy processes that is based on the Wiener-Hopf decomposition. We pursue this idea further by combining their technique with the recently introduced multilevel Monte Carlo methodology. Moreover, we provide here for the first time a theoretical analysis of the new Monte Carlo simulation technique in Kuznetsov et al. (2011) and of its multilevel variant for computing expectations of functions depending on the historical trajectory of a Levy process. We derive rates of convergence for both methods and show that they are uniform with respect to the "jump activity" (e.g. characterised by the Blumenthal-Getoor index). We also present a modified version of the algorithm in Kuznetsov et al. (2011) which combined with the multilevel methodology obtains the optimal rate of convergence for general Levy processes and Lipschitz functionals. This final result is only a theoretical one at present, since it requires independent sampling from a triple of distributions which is currently only possible for a limited number of processes.

0 citations0 reviewsNo code linkedOpen access