We develop a hybrid spatial discretization for the wave equation in second order form, based on high-order accurate finite difference methods and discontinuous Galerkin methods. The hybridization combines computational efficiency of finite difference methods on Cartesian grids and geometrical flexibility of discontinuous Galerkin methods on unstructured meshes. The two spatial discretizations are coupled by a penalty technique at the interface such that the overall semidiscretization satisfies a discrete energy estimate to ensure stability. In addition, optimal convergence is obtained in the sense that when combining a fourth order finite difference method with a discontinuous Galerkin method using third order local polynomials, the overall convergence rate is fourth order. Furthermore, we use a novel approach to derive an error estimate for the semidiscretization by combining the energy method and the normal mode analysis for a corresponding one dimensional model problem. The stability and accuracy analysis are verified in numerical experiments.
Practitioners are often left tuning Metropolis-Hastings algorithms by trial and error or using optimal scaling guidelines to avoid poor empirical performance. We develop lower bounds on the convergence rates of geometrically ergodic accept-reject-based Markov chains (i.e. Metropolis-Hastings, non-reversible Metropolis-Hastings) to study their computational complexity. If the target density concentrates with a parameter $n$ (e.g. Bayesian posterior concentration, Laplace approximations), we show the convergence rate can tend to $1$ exponentially fast if the tuning parameters do not depend carefully on $n$. We show this is the case for random-walk Metropolis in Bayesian logistic regression with Zellner's g-prior when the dimension and sample size $d/n \to \gamma \in (0, 1)$. We focus on more general target densities using a special class of Metropolis-Hastings algorithms with a Gaussian proposal (e.g. random walk and Metropolis-adjusted Langevin algorithms) where we give more general conditions. An application to flat prior Bayesian logistic regression as $n \to \infty$ is studied. We also develop lower bounds in the Wasserstein distances which have become popular in the convergence analysis of high-dimensional MCMC algorithms with similar conclusions.
We introduce an arbitrary order, stabilized finite element method for solving a unique continuation problem subject to the time-harmonic elastic wave equation with variable coefficients. Based on conditional stability estimates we prove convergence rates for the proposed method which take into account the noise level and the polynomial degree. A series of numerical experiments corroborates our theoretical results and explores additional aspects, e.g. how the quality of the reconstruction depends on the geometry of the involved domains. We find that certain convexity properties are crucial to obtain a good recovery of the wave displacement outside the data domain and that higher polynomial orders can be more efficient but also more sensitive to the ill-conditioned nature of the problem.
In this paper, we propose a low-cost, parameter-free, and pressure-robust Stokes solver based on the enriched Galerkin (EG) method with a discontinuous velocity enrichment function. The EG method employs the interior penalty discontinuous Galerkin (IPDG) formulation to weakly impose the continuity of the velocity function. However, the symmetric IPDG formulation, despite of its advantage of symmetry, requires a lot of computational effort to choose an optimal penalty parameter and to compute different trace terms. In order to reduce such effort, we replace the derivatives of the velocity function with its weak derivatives computed by the geometric data of elements. Therefore, our modified EG (mEG) method is a parameter-free numerical scheme which has reduced computational complexity as well as optimal rates of convergence. Moreover, we achieve pressure-robustness for the mEG method by employing a velocity reconstruction operator on the load vector on the right-hand side of the discrete system. The theoretical results are confirmed through numerical experiments with two- and three- dimensional examples.
In this paper we show that conditional graph entropy can be formulated as an alternating minimization problem, which gives rise to a simple iterative algorithm for numerically computing (conditional) graph entropy. The systematic study of alternating minimization problems was initiated by Csisz\'ar and Tusn\'ady. We apply this theory to conditional graph entropy, which was shown to be the minimal rate for a natural functional compression problem with side information at the receiver. This also leads to a new formula which shows that conditional graph entropy is part of a more general framework: the solution of an optimization problem over a convex corner. In the special case of graph entropy (i.e., unconditioned version) this was known due to Csisz\'ar, K\"orner, Lov\'asz, Marton, and Simonyi. In that case the role of the convex corner was played by the so-called vertex packing polytope. In the conditional version it is a more intricate convex body but the function to minimize is the same.
The selection of smoothing parameter is central to the estimation of penalized splines. The best value of the smoothing parameter is often the one that optimizes a smoothness selection criterion, such as generalized cross-validation error (GCV) and restricted likelihood (REML). To correctly identify the global optimum rather than being trapped in an undesired local optimum, grid search is recommended for optimization. Unfortunately, the grid search method requires a pre-specified search interval that contains the unknown global optimum, yet no guideline is available for providing this interval. As a result, practitioners have to find it by trial and error. To overcome such difficulty, we develop novel algorithms to automatically find this interval. Our automatic search interval has four advantages. (i) It specifies a smoothing parameter range where the associated penalized least squares problem is numerically solvable. (ii) It is criterion-independent so that different criteria, such as GCV and REML, can be explored on the same parameter range. (iii) It is sufficiently wide to contain the global optimum of any criterion, so that for example, the global minimum of GCV and the global maximum of REML can both be identified. (iv) It is computationally cheap compared with the grid search itself, carrying no extra computational burden in practice. Our method is ready to use through our recently developed R package gps (>= version 1.1). It may be embedded in more advanced statistical modeling methods that rely on penalized splines.
Simulating propagation of acoustic waves via solving a system of three-coupled first-order linear differential equations using a k-space pseudo-spectral method is popular for biomedical applications, firstly because of availability of an open-source toolbox for implementation of this numerical approach, and secondly because of its efficiency. The k-space pseudo-spectral method is efficient, because it allows coarser computational grids and larger time steps than finite difference and finite element methods for the same accuracy. The goal of this study is to compare this numerical wave solver with an analytical solution to the wave equation using the Green's function for computing propagation of acoustic waves in homogeneous media. This comparison is done in the frequency domain. Using the k-Wave solver, a match to the Green's function is obtained after modifying the approach taken for including mass source in the linearised equation of continuity (conservation of mass) in the associated system of wave equations.
In this work we connect two notions: That of the nonparametric mode of a probability measure, defined by asymptotic small ball probabilities, and that of the Onsager--Machlup functional, a generalized density also defined via asymptotic small ball probabilities. We show that in a separable Hilbert space setting and under mild conditions on the likelihood, the modes of a Bayesian posterior distribution based upon a Gaussian prior agree with the minimizers of its Onsager--Machlup functional. We apply this result to inverse problems and derive conditions on the forward mapping under which this variational characterization of posterior modes holds. Our results show rigorously that in the limit case of infinite-dimensional data corrupted by additive Gaussian or Laplacian noise, nonparametric MAP estimation is equivalent to Tikhonov--Phillips regularization. In comparison with the work of Dashti, Law, Stuart, and Voss (2013), the assumptions on the likelihood are relaxed so that they cover in particular the important case of Gaussian process noise. We illustrate our results by applying them to a severely ill-posed linear problem with Laplacian noise, where we express the MAP estimator analytically and study its rate of convergence.
The closure principle is fundamental in multiple testing and has been used to derive many efficient procedures with familywise error rate control. However, it is often not suitable for modern research, as more flexible multiple testing settings are considered where not all hypotheses are known at the beginning of the evaluation. In this paper, we focus on online multiple testing where a possibly infinite sequence of hypotheses is tested over time. At each step, it must be decided on the current hypothesis without having any information about the hypotheses that have not been tested yet. Our main contribution is a new online closure principle which ensures that the resulting closed procedure can be applied in the online setting. We prove that any familywise error rate (FWER) controlling online procedure can be derived by this online closure principle. In addition, we demonstrate how short-cuts of these online closed procedures can be obtained under a suitable consonance property and apply the results in order to construct new online multiple testing methods. Finally, the new online closure principle is used to derive an improvement of the currently most promising online procedure with FWER control, the ADDIS-Spending under local dependence.
In two and three dimension we analyze discontinuous Galerkin methods for the acoustic problem. The acoustic fluid that we consider on this paper is inviscid, leading to a linear eigenvalue problem. The acoustic problem is written, in first place, in terms of the displacement. Under the approach of the non-compact operators theory, we prove convergence and error estimates for the method when the displacement formulation is considered. We analyze the influence of the stabilization parameter on the computation of the spectrum, where spurious eigenmodes arise when this parameter is not correctly chosen. Alternatively we present the formulation depending only on the pressure, comparing the performance of the DG methods with the pure displacement formulation. Computationally, we study the influence of the stabilization parameter on the arising of spurious eigenvalues when the spectrum is computed. Also, we report tests in two and three dimensions where convergence rates are reported, together with a comparison between the displacement and pressure formulations for the proposed DG methods.
Classic algorithms and machine learning systems like neural networks are both abundant in everyday life. While classic computer science algorithms are suitable for precise execution of exactly defined tasks such as finding the shortest path in a large graph, neural networks allow learning from data to predict the most likely answer in more complex tasks such as image classification, which cannot be reduced to an exact algorithm. To get the best of both worlds, this thesis explores combining both concepts leading to more robust, better performing, more interpretable, more computationally efficient, and more data efficient architectures. The thesis formalizes the idea of algorithmic supervision, which allows a neural network to learn from or in conjunction with an algorithm. When integrating an algorithm into a neural architecture, it is important that the algorithm is differentiable such that the architecture can be trained end-to-end and gradients can be propagated back through the algorithm in a meaningful way. To make algorithms differentiable, this thesis proposes a general method for continuously relaxing algorithms by perturbing variables and approximating the expectation value in closed form, i.e., without sampling. In addition, this thesis proposes differentiable algorithms, such as differentiable sorting networks, differentiable renderers, and differentiable logic gate networks. Finally, this thesis presents alternative training strategies for learning with algorithms.