Divakar Viswanath

Divakar Viswanath

Spectral integration was deployed by Orszag and co-workers (1977, 1980, 1981) to obtain stable and efficient solvers for the incompressible Navier-Stokes equation in rectangular geometries. Two methods in current use for channel flow and plane Couette flow, namely, Kleiser-Schumann (1980) and Kim-Moin-Moser (1977), rely on the same technique. In its current form, the technique of spectral integration, as applied to the Navier-Stokes equations, is dominated by rounding errors at higher Reynolds numbers which would otherwise be within reach. In this article, we derive a number of versions of spectral integration and explicate their properties, with a view to extending the Kleiser-Schumann and Kim-Moin-Moser algorithms to higher Reynolds numbers. More specifically, we show how spectral integration matrices that are banded, but bordered by dense rows, can be reduced to purely banded matrices. Key properties, such as the accuracy of spectral integration even when Green's functions are not resolved by the underlying grid, the accuracy of spectral integration in spite of ill-conditioning of underlying linear systems, and the accuracy of derivatives, are thoroughly explained.

Burhan Sadiq, Divakar Viswanath

Let $z_{1},\ldots,z_{K}$ be distinct grid points. If $f_{k,0}$ is the prescribed value of a function at the grid point $z_{k}$, and $f_{k,r}$ the prescribed value of the $r$\foreignlanguage{american}{-th} derivative, for $1\leq r\leq n_{k}-1$, the Hermite interpolant is the unique polynomial of degree $N-1$ ($N=n_{1}+\cdots+n_{K}$) which interpolates the prescribed function values and function derivatives. We obtain another derivation of a method for Hermite interpolation recently proposed by Butcher et al. {[}\emph{Numerical Algorithms, vol. 56 (2011), p. 319-347}{]}. One advantage of our derivation is that it leads to an efficient method for updating the barycentric weights. If an additional derivative is prescribed at one of the interpolation points, we show how to update the barycentric coefficients using only $\mathcal{O}\left(N\right)$ operations. Even in the context of confluent Newton series, a comparably efficient and general method to update the coefficients appears not to be known. If the method is properly implemented, it computes the barycentric weights with fewer operations than other methods and has very good numerical stability even when derivatives of high order are involved. We give a partial explanation of its numerical stability.

Burhan Sadiq, Divakar Viswanath

Let $z_{1},z_{2},...,z_{N}$ be a sequence of distinct grid points. A finite difference formula approximates the $m$-th derivative $f^{(m)}(0)$ as $\sum w_{k}f(z_{k})$, with $w_{k}$ being the weights. We derive an algorithm for finding the weights $w_{k}$ which is an improvement of an algorithm of Fornberg (\emph{Mathematics of Computation}, vol. 51 (1988), p. 699-706). This algorithm uses fewer arithmetic operations than that of Fornberg by a factor of $4/(5m+5)$ while being equally accurate. The algorithm that we derive computes finite difference weights accurately even when $m$, the order of the derivative, is as high as 16. In addition, the algorithm generalizes easily to the efficient computation of spectral differentiation matrices. The order of accuracy of the finite difference formula for $f^{(m)}(0)$ with grid points $hz_{k}$, $1\leq k\leq N$, is typically $\mathcal{O}(h^{N-m})$. However, the most commonly used finite difference formulas have an order of accuracy that is higher than the typical. For instance, the centered difference approximation $(f(h)-2f(0)+f(-h))/h^{2}$ to $f"(0)$ has an order of accuracy equal to 2 not 1. Even unsymmetric finite difference formulas can exhibit such superconvergence or boosted order of accuracy, as shown by the explicit algebraic condition that we derive. If the grid points are real, we prove a basic result stating that the order of accuracy can never be boosted by more than 1.

Raymundo Navarrete, Divakar Viswanath

This article considers the numerical inversion of the power series $p(x)=1+b_{1}x+b_{2}x^{2}+\cdots$ to compute the inverse series $q(x)$ satisfying $p(x)q(x)=1$. Numerical inversion is a special case of triangular back-substitution, which has been known for its beguiling numerical stability since the classic work of Wilkinson (1961). We prove the numerical stability of inversion of power series and obtain bounds on numerical error. A range of examples show these bounds to be quite good. When $p(x)$ is a polynomial and $x=a$ is a root with $p(a)=0$, we show that root deflation via the simple division $p(x)/(x-a)$ can trigger instabilities relevant to polynomial root finding and computation of finite-difference weights. When $p(x)$ is a polynomial, the accuracy of the computed inverse $q(x)$ is connected to the pseudozeros of $p(x)$.

Divakar Viswanath

The Chebyshev points are commonly used for spectral differentiation in non-periodic domains. The rounding error in the Chebyshev approximation to the $n$-the derivative increases at a rate greater than $n^{2m}$ for the $m$-th derivative. The mapping technique of Kosloff and Tal-Ezer (\emph{J. Comp. Physics}, vol. 104 (1993), p. 457-469) ameliorates this increase in rounding error. We show that the argument used to justify the choice of the mapping parameter is substantially incomplete. We analyze rounding error as well as discretization error and give a more complete argument for the choice of the mapping parameter. If the discrete cosine transform is used to compute derivatives, we show that a different choice of the mapping parameter yields greater accuracy.

Raymundo Navarrete, Divakar Viswanath

Prediction of dynamical time series with additive noise using support vector machines or kernel based regression has been proved to be consistent for certain classes of discrete dynamical systems. Consistency implies that these methods are effective at computing the expected value of a point at a future time given the present coordinates. However, the present coordinates themselves are noisy, and therefore, these methods are not necessarily effective at removing noise. In this article, we consider denoising and prediction as separate problems for flows, as opposed to discrete time dynamical systems, and show that the use of smooth splines is more effective at removing noise. Combination of smooth splines and kernel based regression yields predictors that are more accurate on benchmarks typically by a factor of 2 or more. We prove that kernel based regression in combination with smooth splines converges to the exact predictor for time series extracted from any compact invariant set of any sufficiently smooth flow. As a consequence of convergence, one can find examples where the combination of kernel based regression with smooth splines is superior by even a factor of $100$. The predictors that we compute operate on delay coordinate data and not the full state vector, which is typically not observable.