Search

We consider problems of minimizing functionals $\mathcal{F}$ of probability measures on the Euclidean space. To propose an accelerated gradient descent algorithm for such problems, we consider gradient flow of transport maps that give push-forward measures of an initial measure. Then we propose a deterministic accelerated algorithm by extending Nesterov's acceleration technique with momentum. This algorithm do not based on the Wasserstein geometry. Furthermore, to estimate the convergence rate of the accelerated algorithm, we introduce new convexity and smoothness for $\mathcal{F}$ based on transport maps. As a result, we can show that the accelerated algorithm converges faster than a normal gradient descent algorithm. Numerical experiments support this theoretical result., Comment: 31 pages

Although numerous studies have focused on normal Besov spaces, limited studies have been conducted on exponentially weighted Besov spaces. Therefore, we define exponentially weighted Besov space $VB_{p,q}^{\delta,w}(\mathbb{R}^d)$ whose smoothness includes normal Besov spaces, Besov spaces with dominating mixed smoothness, and their interpolation. Furthermore, we obtain wavelet characterization of $VB_{p,q}^{\delta,w}(\mathbb{R}^d)$. Next, approximation formulas such as sparse grids are derived using the determined formula. The results of this study are expected to provide considerable insight into the application of exponentially weighted Besov spaces with mixed smoothness., Comment: 35 pages, 2 figures

In this paper, we claim that a common underlying structure--a skeleton structure--is present behind discrete-time quantum walks (QWs) on a one-dimensional lattice with a homogeneous coin matrix. This skeleton structure is independent of the initial state, and partially, even of the coin matrix. This structure is best interpreted in the context of quantum-walk-replicating random walks (QWRWs), i.e., random walks that replicate the probability distribution of quantum walks, where this newly found structure acts as a simplified formula for the transition probability. Additionally, we construct a random walk whose transition probabilities are defined by the skeleton structure and demonstrate that the resultant properties of the walkers are similar to both the original QWs and QWRWs., Comment: 26 pages, 9 figures

We consider a class of particle systems which appear in various applications such as approximation theory, plasticity, potential theory and space-filling designs. The positions of the particles on the real line are described as a global minimum of an interaction energy, which consists of a nonlocal, repulsive interaction part and a confining part. Motivated by the applications, we cover non-standard scenarios in which the confining potential weakens as the number of particles increases. This results in a large area over which the particles spread out. Our aim is to approximate the particle interaction energy by a corresponding continuum interacting energy. Our main results are bounds on the corresponding energy difference and on the difference between the related potential values. We demonstrate that these bounds are useful to problems in approximation theory and plasticity. The proof of these bounds relies on convexity assumptions on the interaction and confining potentials. It combines recent advances in the literature with a new upper bound on the minimizer of the continuum interaction energy., Comment: 28 pages

Based on the Sinc approximation combined with the tanh transformation, Haber derived an approximation formula for numerical indefinite integration over the finite interval (-1, 1). The formula uses a special function for the basis functions. In contrast, Stenger derived another formula, which does not use any special function but does include a double sum. Subsequently, Muhammad and Mori proposed a formula, which replaces the tanh transformation with the double-exponential transformation in Haber's formula. Almost simultaneously, Tanaka et al. proposed another formula, which was based on the same replacement in Stenger's formula. As they reported, the replacement drastically improves the convergence rate of Haber's and Stenger's formula. In addition to the formulas above, Stenger derived yet another indefinite integration formula based on the Sinc approximation combined with the tanh transformation, which has an elegant matrix-vector form. In this paper, we propose the replacement of the tanh transformation with the double-exponential transformation in Stenger's second formula. We provide a theoretical analysis as well as a numerical comparison.

The Pairwise Conditional Gradients (PCG) algorithm is a powerful extension of the Frank-Wolfe algorithm leading to particularly sparse solutions, which makes PCG very appealing for problems such as sparse signal recovery, sparse regression, and kernel herding. Unfortunately, PCG exhibits so-called swap steps that might not provide sufficient primal progress. The number of these bad steps is bounded by a function in the dimension and as such known guarantees do not generalize to the infinite-dimensional case, which would be needed for kernel herding. We propose a new variant of PCG, the so-called Blended Pairwise Conditional Gradients (BPCG). This new algorithm does not exhibit any swap steps, is very easy to implement, and does not require any internal gradient alignment procedures. The convergence rate of BPCG is basically that of PCG if no drop steps would occur and as such is no worse than PCG but improves and provides new rates in many cases. Moreover, we observe in the numerical experiments that BPCG's solutions are much sparser than those of PCG. We apply BPCG to the kernel herding setting, where we derive nice quadrature rules and provide numerical results demonstrating the performance of our method.

The kernel herding algorithm is used to construct quadrature rules in a reproducing kernel Hilbert space (RKHS). While the computational efficiency of the algorithm and stability of the output quadrature formulas are advantages of this method, the convergence speed of the integration error for a given number of nodes is slow compared to that of other quadrature methods. In this paper, we propose a modified kernel herding algorithm whose framework was introduced in a previous study and aim to obtain sparser solutions while preserving the advantages of standard kernel herding. In the proposed algorithm, the negative gradient is approximated by several vertex directions, and the current solution is updated by moving in the approximate descent direction in each iteration. We show that the convergence speed of the integration error is directly determined by the cosine of the angle between the negative gradient and approximate gradient. Based on this, we propose new gradient approximation algorithms and analyze them theoretically, including through convergence analysis. In numerical experiments, we confirm the effectiveness of the proposed algorithms in terms of sparsity of nodes and computational efficiency. Moreover, we provide a new theoretical analysis of the kernel quadrature rules with fully-corrective weights, which realizes faster convergence speeds than those of previous studies.

We propose a method for generating nodes for kernel quadrature by a point-wise gradient descent method. For kernel quadrature, most methods for generating nodes are based on the worst case error of a quadrature formula in a reproducing kernel Hilbert space corresponding to the kernel. In typical ones among those methods, a new node is chosen among a candidate set of points in each step by an optimization problem with respect to a new node. Although such sequential methods are appropriate for adaptive quadrature, it is difficult to apply standard routines for mathematical optimization to the problem. In this paper, we propose a method that updates a set of points one by one with a simple gradient descent method. To this end, we provide an upper bound of the worst case error by using the fundamental solution of the Laplacian on $\mathbf{R}^{d}$. We observe the good performance of the proposed method by numerical experiments., Comment: 21 pages, 12 figures

In this paper, we investigate application of mathematical optimization to construction of a cubature formula on Wiener space, which is a weak approximation method of stochastic differential equations introduced by Lyons and Victoir (Cubature on Wiener Space, Proc. R. Soc. Lond. A 460, 169--198). After giving a brief review of the cubature theory on Wiener space, we show that a cubature formula of general dimension and degree can be obtained through a Monte Carlo sampling and linear programming. This paper also includes an extension of stochastic Tchakaloff's theorem, which technically yields the proof of our primary result., Comment: 25 pages

We construct approximate Fekete point sets for kernel-based interpolation by maximising the determinant of a kernel Gram matrix obtained via truncation of an orthonormal expansion of the kernel. Uniform error estimates are proved for kernel interpolants at the resulting points. If the kernel is Gaussian we show that the approximate Fekete points in one dimension are the solution to a convex optimisation problem and that the interpolants converge with a super-exponential rate. Numerical examples are provided for the Gaussian kernel., Comment: Accepted for publication in Numerical Algorithms

We investigate the approximation formulas that were proposed by Tanaka & Sugihara (2019), in weighted Hardy spaces, which are analytic function spaces with certain asymptotic decay. Under the criterion of minimum worst error of $n$-point approximation formulas, we demonstrate that the formulas are nearly optimal. We also obtain the upper bounds of the approximation errors that coincide with the existing heuristic bounds in asymptotic order by duality theorem for the minimization problem of potential energy., Comment: 17 pages

The double exponential formula, or the DE formula, is a high-precision integration formula using a change of variables called a DE transformation; whereas there is a disadvantage that it is sensitive to singularities of an integrand near the real axis. To overcome this disadvantage, Slevinsky and Olver (SIAM J. Sci. Comput., 2015) attempted to improve it by constructing conformal maps based on the locations of singularities. Based on their ideas, we construct a new transformation formula. Our method employs special types of the Schwarz-Christoffel transformations for which we can derive their explicit form. Then, the new transformation formula can be regarded as a generalization of the DE transformations. We confirm its effectiveness by numerical experiments.

We propose algorithms to take point sets for kernel-based interpolation of functions in reproducing kernel Hilbert spaces (RKHSs) by convex optimization. We consider the case of kernels with the Mercer expansion and propose an algorithm by deriving a second-order cone programming (SOCP) problem that yields $n$ points at one sitting for a given integer $n$. In addition, by modifying the SOCP problem slightly, we propose another sequential algorithm that adds an arbitrary number of new points in each step. Numerical experiments show that in several cases the proposed algorithms compete with the $P$-greedy algorithm, which is known to provide nearly optimal points., Comment: 31 pages. The programs for the numerical computation in this article are available on https://github.com/KeTanakaN/mat_points_interp_rkhs

In this paper, we address the problem of counting integer points in a rational polytope described by $P(y) = \{ x \in \mathbb{R}^m \colon Ax = y, x \geq 0\}$, where $A$ is an $n \times m$ integer matrix and $y$ is an $n$-dimensional integer vector. We study the Z-transformation approach initiated by Brion-Vergne, Beck, and Lasserre-Zeron from the numerical analysis point of view, and obtain a new algorithm on this problem: If $A$ is nonnegative, then the number of integer points in $P(y)$ can be computed in $O(\mathrm{poly} (n,m, \|y\|_\infty) (\|y\|_\infty + 1)^n)$ time and $O(\mathrm{poly} (n,m, \|y\|_\infty))$ space.This improves, in terms of space complexity, a naive DP algorithm with $O((\|y\|_\infty + 1)^n)$-size DP table. Our result is based on the standard error analysis to the numerical contour integration for the inverse Z-transform, and establish a new type of an inclusion-exclusion formula for integer points in $P(y)$. We apply our result to hypergraph $b$-matching, and obtain a $O(\mathrm{poly}( n,m,\|b\|_\infty) (\|b\|_\infty +1)^{(1-1/k)n})$ time algorithm for counting $b$-matchings in a $k$-partite hypergraph with $n$ vertices and $m$ hyperedges. This result is viewed as a $b$-matching generalization of the classical result by Ryser for $k=2$ and its multipartite extension by Bj{\"o}rklund-Husfeldt., Comment: 13 pages

We propose a simple and effective method for designing approximation formulas for weighted analytic functions. We consider spaces of such functions according to weight functions expressing the decay properties of the functions. Then, we adopt the minimum worst error of the $n$-point approximation formulas in each space for characterizing the optimal sampling points for the approximation. In order to obtain approximately optimal sampling points, we consider minimization of a discrete energy related to the minimum worst error. Consequently, we obtain an approximation formula and its theoretical error estimate in each space. In addition, from some numerical experiments, we observe that the formula generated by the proposed method outperforms the corresponding formula derived with sinc approximation, which is near optimal in each space., Comment: 27 pages, 15 figures. The programs for the numerical computations in this article are available on https://github.com/KeTanakaN/mat_disc_ener_opt_approx

We propose an optimal approximation formula for analytic functions that are defined on a complex region containing the real interval $(-1,1)$ and possibly have algebraic singularities at the endpoints of the interval. As a space of such functions,we consider a Hardy space with the weight given by $w_{\mu}(z) = (1-z^{2})^{\mu/2}$ for $\mu > 0$, and formulate the optimality of an approximation formula for the functions in the space. Then, we propose an optimal approximation formula for the space for any $\mu > 0$ as opposed to existing results with the restriction $0 < \mu < \mu_{\ast}$ for a certain constant $\mu_{\ast}$. We also provide the results of numerical experiments to show the performance of the proposed formula., Comment: 24 pages, 7 figures

We propose a method for designing accurate interpolation formulas on the real axis for the purpose of function approximation in weighted Hardy spaces. In particular, we consider the Hardy space of functions that are analytic in a strip region around the real axis, being characterized by a weight function $w$ that determines the decay rate of its elements in the neighborhood of infinity. Such a space is considered as a set of functions that are transformed by variable transformations that realize a certain decay rate at infinity. Popular examples of such transformations are given by the single exponential (SE) and double exponential (DE) transformations for the SE-Sinc and DE-Sinc formulas, which are very accurate owing to the accuracy of sinc interpolation in the weighted Hardy spaces with single and double exponential weights $w$, respectively. However, it is not guaranteed that the sinc formulas are optimal in weighted Hardy spaces, although Sugihara has demonstrated that they are near optimal. An explicit form for an optimal approximation formula has only been given in weighted Hardy spaces with SE weights of a certain type. In general cases, explicit forms for optimal formulas have not been provided so far. We adopt a potential theoretic approach to obtain almost optimal formulas in weighted Hardy spaces in the case of general weight functions $w$. We formulate the problem of designing an optimal formula in each space as an optimization problem written in terms of a Green potential with an external field. By solving the optimization problem numerically, we obtain an almost optimal formula in each space. Furthermore, some numerical results demonstrate the validity of this method. In particular, for the case of a DE weight, the formula designed by our method outperforms the DE-Sinc formula., Comment: 40 pages, 15 figures, Matlab codes for this article are available on https://github.com/KeTanakaN/mat_pot_approx

In this paper, we propose a fast and accurate numerical method based on Fourier transform to solve Kolmogorov forward equations of symmetric scalar L\'evy processes. The method is based on the accurate numerical formulas for Fourier transform proposed by Ooura. These formulas are combined with nonuniform fast Fourier transform (FFT) and fractional FFT to speed up the numerical computations. Moreover, we propose a formula for numerical indefinite integration on equispaced grids as a component of the method. The proposed integration formula is based on the sinc-Gauss sampling formula, which is a function approximation formula. This integration formula is also combined with the FFT. Therefore, all steps of the proposed method are executed using the FFT and its variants. The proposed method allows us to be free from some special treatments for a non-smooth initial condition and numerical time integration. The numerical solutions obtained by the proposed method appeared to be exponentially convergent on the interval if the corresponding exact solutions do not have sharp cusps. Furthermore, the real computational times are approximately consistent with the theoretical estimates., Comment: 24 pages, 15 figures, Matlab codes for this article are available on https://github.com/KeTanakaN/mat_symLevy_FT_codes

The maximum entropy principle is a powerful tool for solving underdetermined inverse problems. This paper considers the problem of discretizing a continuous distribution, which arises in various applied fields. We obtain the approximating distribution by minimizing the Kullback-Leibler information (relative entropy) of the unknown discrete distribution relative to an initial discretization based on a quadrature formula subject to some moment constraints. We study the theoretical error bound and the convergence of this approximation method as the number of discrete points increases. We prove that (i) the theoretical error bound of the approximate expectation of any bounded continuous function has at most the same order as the quadrature formula we start with, and (ii) the approximate discrete distribution weakly converges to the given continuous distribution. Moreover, we present some numerical examples that show the advantage of the method and apply to numerically solving an optimal portfolio problem., Comment: 20 pages, 14 figures

In this paper, we consider a method for fast numerical computation of the Fourier transform of a slowly decaying function with given accuracy in given ranges of the frequency. In these decades, some useful formulas for the Fourier transform are proposed to recover difficulty of the computation due to the slow decay and the oscillation of the integrand. In particular, Ooura proposed formulas with continuous Euler transformation and showed their effectiveness. It is, however, also reported that errors of them become large outside some ranges of the frequency. Then, for an illustrating representative of the formulas, we choose parameters in the formula based on its error analysis to compute the Fourier transform with given accuracy in given ranges of the frequency. Furthermore, combining the formula and fractional FFT, a generalization of the fast Fourier transform (FFT), we execute the computation in the same order of computation time as the FFT., Comment: 15 pages, 13 figures

We present the convergence rates and the explicit error bounds of Hill's method, which is a numerical method for computing the spectra of ordinary differential operators with periodic coefficients. This method approximates the operator by a finite dimensional matrix. On the assumption that the operator is selfadjoint, it is shown that, under some conditions, we can obtain the convergence rates of eigenvalues with respect to the dimension and the explicit error bounds. Numerical examples demonstrate that we can verify these conditions using Gershgorin's theorem for some real problems. Main theorems are proved using the Dunford integrals which project an eigenvector to the corresponding eigenspace., Comment: 22 pages, 4 figures