Fast and accurate tensor approximation of a multivariate convolution with linear scaling in dimension

In the present paper we present the tensor-product approximation of a multidimensional convolution transform discretized via a collocation-projection scheme on uniform or composite refined grids. Examples of convolving kernels are provided by the classical Newton, Slater (exponential) and Yukawa potentials, 1/‖x‖, and with x∈R^d. For piecewise constant elements on the uniform grid of size n^d, we prove quadratic convergence O(h²) in the mesh parameter h=1/n, and then justify the Richardson extrapolation method on a sequence of grids that improves the order of approximation up to O(h³). A fast algorithm of complexity O(dR₁R₂nlogn) is described for tensor-product convolution on uniform/composite grids of size n^d, where R₁,R₂ are tensor ranks of convolving functions. We also present the tensor-product convolution scheme in the two-level Tucker canonical format and discuss the consequent rank reduction strategy. Finally, we give numerical illustrations confirming: (a) the approximation theory for convolution schemes of order O(h²) and O(h³); (b) linear-logarithmic scaling of 1D discrete convolution on composite grids; (c) linear-logarithmic scaling in n of our tensor-product convolution method on an n×n×n grid in the range n≤16384.     

Sign regular matrices and Neville elimination

A pivoting strategy of O(n) operations for the Neville elimination of n × n nonsingular sign regular matrices is introduced. Among other nice properties, it is proved that it preserves sign regularity. It is also shown its relationship with scaled partial pivoting strategies for Neville elimination.     

A sparse counterpart of Reichel and Gragg’s package QRUP

We describe how to maintain the triangular factor of a sparse QR factorization when columns are added and deleted and Q cannot be stored for sparsity reasons. The updating procedures could be thought of as a sparse counterpart of Reichel and Gragg’s package QRUP. They allow us to solve a sequence of sparse linear least squares subproblems in which each matrix B_k is an independent subset of the columns of a fixed matrix A, and B_k+1 is obtained by adding or deleting one column. Like Coleman and Hulbert [T. Coleman, L. Hulbert, A direct active set algorithm for large sparse quadratic programs with simple bounds, Math. Program. 45 (1989) 373-406], we adapt the sparse direct methodology of Björck and Oreborn of the late 1980s, but without forming A^TA, which may be only positive semidefinite. Our Matlab 5 implementation works with a suitable row and column numbering within a static triangular sparsity pattern that is computed in advance by symbolic factorization of A^TA and preserved with placeholders.     

Covering times of random walks on bounded degree trees and other graphs

The motivating problem for this paper is to find the expected covering time of a random walk on a balanced binary tree withn vertices. Previous upper bounds for general graphs ofO(|V| |E|)⁽¹⁾ andO(|V| |E|/d _min)⁽²⁾ imply an upper bound ofO(n ²). We show an upper bound on general graphs ofO( |E| log |V|), which implies an upper bound ofO(n log² n). The previous lower bound was (|V| log |V|) for trees.⁽²⁾ In our main result, we show a lower bound of (|V| (log _{d max} |V|)²) for trees, which yields a lower bound of (n log² n). We also extend our techniques to show an upper bound for general graphs ofO(max{E T_i} log |V|).     

NA随机场对数律的收敛速度

设d是一个正整数, N ^d是d -维正整数格点.设{X_n , n∈N ^d} 是一同分布的负相伴随机场, 记S_n =∑_{k≤ n} X_k, S_n^(k)=S_n-X_k, 如果r >2, EX₁ = 0 和σ²= Var(X₁}, 则存在一个正数M:=100√(r-2)(1+σ²)使得下列条件等价 (I)E |X₁|^r (log|X₁|)^d-1-r/2 <∞; (II)∑_{n∈ N^d} |n|^r/2-2P(max_{1≤ k≤ n} |S_n^(k)|≥ (2^d+1 )ε√|n| log |n |) <∞,∨ε > M; (III)∑_{n∈N ^d} |n|^r/2-2P(max_{1≤ k≤n} |S_k |≥ε√| n} log| n |) <∞,∨ε > M. (III)\ \ $\sum\limits_{{{\bf n}}\in {{\cal N}}^{d}} |n|^{r/2-2} P(\max\limits_{{\bf 1}\leq{\bf k}\leq{\bf n}}|S_{{\bf k}}|\geq \varepsilon \sqrt{|{\bf n}|\log |{\bf n}|})<\infty$, $\forall\varepsilon>M$.     

A fast and simple randomized parallel algorithm for the maximal independent set problem

A simple parallel randomized algorithm to find a maximal independent set in a graph G = (V, E) on n vertices is presented. Its expected running time on a concurrent-read concurrent-write PRAM with O(|E|d_max) processors is O(log n), where d_max denotes the maximum degree. On an exclusive-read exclusive-write PRAM with O(|E|) processors the algorithm runs in O(log²n). Previously, an O(log⁴n) deterministic algorithm was given by Karp and Wigderson for the EREW-PRAM model. This was recently (independently of our work) improved to O(log²n) by M. Luby. In both cases randomized algorithms depending on pairwise independent choices were turned into deterministic algorithms. We comment on how randomized combinatorial algorithms whose analysis only depends on d-wise rather than fully independent random choices (for some constant d) can be converted into deterministic algorithms. We apply a technique due to A. Joffe (1974) and obtain deterministic construction in fast parallel time of various combinatorial objects whose existence follows from probabilistic arguments.     

The analysis of a nested dissection algorithm

Summary Nested dissection is an algorithm invented by Alan George for preserving sparsity in Gaussian elimination on symmetric positive definite matrices. Nested dissection can be viewed as a recursive divide-and-conquer algorithm on an undirected graph; it usesseparators in the graph, which are small sets of vertices whose removal divides the graph approximately in half. George and Liu gave an implementation of nested dissection that used a heuristic to find separators. Lipton and Tarjan gave an algorithm to findn ^1/2-separators in planar graphs and two-dimensional finite element graphs, and Lipton, Rose, and Tarjan used these separators in a modified version of nested dissection, guaranteeing bounds ofO (n logn) on fill andO(n ^3/2) on operation count. We analyze the combination of the original George-Liu nested dissection algorithm and the Lipton-Tarjan planar separator algorithm. This combination is interesting because it is easier to implement than the Lipton-Rose-Tarjan version, especially in the framework of existïng sparse matrix software. Using some topological graph theory, we proveO(n logn) fill andO(n ^3/2) operation count bounds for planar graphs, twodimensional finite element graphs, graphs of bounded genus, and graphs of bounded degree withn ^1/2-separators. For planar and finite element graphs, the leading constant factor is smaller than that in the Lipton-Rose-Tarjan analysis. We also construct a class of graphs withn ^1/2-separators for which our algorithm does not achieve anO(n logn) bound on fill.The work of this author was supported in part by the Hertz Foundation under a graduate fellowship and by the National Science Foundation under Grant MCS 82-02948The work of this author was supported in part by the National Science Foundation under Grant MCS 78-26858 and by the Office of Naval Research under Contract N00014-76-C-0688     

Fast inversion of vandermonde-like matrices involving orthogonal polynomials

Let {q} _j ^=0n–1 be a family of polynomials that satisfy a three-term recurrence relation and let {t _k } _k ⁼¹ⁿ be a set of distinct nodes. Define the Vandermonde-like matrixW _n =[w _jk ] _k,j ⁼¹ⁿ ,w _jk =q _j–1(t _k ). We describe a fast algorithm for computing the elements of the inverse ofW _n inO(n ²) arithmetic operations. Our algorithm generalizes a scheme presented by Traub [22] for fast inversion of Vandermonde matrices. Numerical examples show that our scheme often yields higher accuracy than the LINPACK subroutine SGEDI for inverting a general matrix. SGEDI uses Gaussian elimination with partial pivoting and requiresO(n ³) arithmetic operations.Dedicated to Gene H. Golub on his 60th birthdayResearch supported by NSF grant DMS-9002884.     

Factorizing symmetric indefinite matrices

The LDL^T factorization of a symmetric indefinite matrix, although efficient computationally, may not exist and can be unstable in the presence of round off error. The use of block diagonal 2×2 pivots is attractive, but there are some difficulties in determining an efficient and stable pivot strategy. Previous suggestions have required O(n>³) operations (either multiplications or comparisons) just to implement the pivot strategy. A new strategy is described which in practice only requires O(n²) operations. Indeed, the effort required by this pivot strategy is less than that required when using partial pivoting with an unsymmetric LU factorization, which is the usual way of factorizing indefinite matrices.     

An elementary operator with log-hyponormal, p-hyponormal entries

A Hilbert space operator A∈B(H) is p-hyponormal, A∈(p-H), if |A^∗|^2p?|A|^2p; an invertible operator A∈B(H) is log-hyponormal, A∈(?-H), if log(TT^∗)?log(T^∗T). Let d_AB=δ_AB or ?_AB, where δ_AB∈B(B(H)) is the generalised derivation δ_AB(X)=AX-XB and ?_AB∈B(B(H)) is the elementary operator ?_AB(X)=AXB-X. It is proved that if A,B^∗∈(?-H)∪(p-H), then, for all complex λ, , the ascent of (d_AB-λ)?1, and d_AB satisfies the range-kernel orthogonality inequality ‖X‖?‖X-(d_AB-λ)Y‖ for all X∈(d_AB-λ)^-1(0) and Y∈B(H). Furthermore, isolated points of σ(d_AB) are simple poles of the resolvent of d_AB. A version of the elementary operator E(X)=A₁XA₂-B₁XB₂ and perturbations of d_AB by quasi-nilpotent operators are considered, and Weyl’s theorem is proved for d_AB.     

Comparison of the best uniform approximations of analytic functions in the disk and on its boundary

Denote by C _A the set of functions that are analytic in the disk |z| < 1 and continuous on its closure |z| ≤ 1; let ? _n , n = 0, 1, 2, ..., be the set of rational functions of degree at most n. Denote by R _n (f) (R _n (f) _A ) the best uniform approximation of a function f ∈ C _A on the circle |z| = 1 (in the disk |z| ≤ 1) by the set ? _n . The following equality is proved for any n ≥ 1: sup{R _n (f) _A /R _n (f): f ∈ C _A ? ? _n } = 2. We also consider a similar problem of comparing the best approximations of functions in C _A by polynomials and trigonometric polynomials.     

Perturbation bounds for the Cholesky andQR factorizations

LetA, A+E be Hermitian positive definite matrices. Suppose thatA=LL ^H andA+E=(L+G)(L+G)^H are the Cholesky factorizations ofA andA+E, respectively. In this paper lower bounds and upper bounds on |G|/|L| in terms of |E|/|A| are given. Moreover, perturbation bounds are given for the QR factorization of a complexm ×n matrixA of rankn.This research was supported by the National Science Foundation of China and the Department of Mathematics of Linköping University in Sweden.     

A fast implicit QR eigenvalue algorithm for companion matrices

An implicit version of the shifted QR eigenvalue algorithm given in Bini et al. [D.A. Bini, Y. Eidelman, I. Gohberg, L. Gemignani, SIAM J. Matrix Anal. Appl. 29(2) (2007) 566-585] is presented for computing the eigenvalues of an n×n companion matrix using O(n²) flops and O(n) memory storage. Numerical experiments and comparisons confirm the effectiveness and the stability of the proposed method.     

On operator norms of submatrices

Both of the following conditions are equivalent to the absoluteness of a norm ν in $\text{C}$ⁿ: (1) for all n×n diagonal matrices D=(d_k), the subordinate operator norm N_ν(D)=max_k|d_k|; (2) for all n×n matrices A, N_ν(A) ?N_ν(|A|). These conditions are modified for partitioned matrices by replacing absolute values with norms of blocks. A generalization of absoluteness is thus obtained.     

Some lower bounds for the spectral radius of matrices using traces

Let A be an n×n matrix with eigenvalues λ₁,λ₂,…,λ_n, and let m be an integer satisfying rank(A)?m?n. If A is real, the best possible lower bound for its spectral radius in terms of m, trA and trA² is obtained. If A is any complex matrix, two lower bounds for are compared, and furthermore a new lower bound for the spectral radius is given only in terms of trA,trA²,‖A‖,‖A^∗A-AA^∗‖,n and m.     

Defect indices of powers of a contraction

Let A be a contraction on a Hilbert space H. The defect index d_A of A is, by definition, the dimension of the closure of the range of I-A^∗A. We prove that (1) d_Aⁿ?nd_A for all n?0, (2) if, in addition, Aⁿ converges to 0 in the strong operator topology and d_A=1, then d_Aⁿ=n for all finite n,0?n?dimH, and (3) d_A=d_A^∗ implies d_Aⁿ=d_A^n∗ for all n?0. The norm-one index k_A of A is defined as sup{n?0:‖Aⁿ‖=1}. When dimH=m<∞, a lower bound for k_A was obtained before: k_A?(m/d_A)-1. We show that the equality holds if and only if either A is unitary or the eigenvalues of A are all in the open unit disc, d_A divides m and d_Aⁿ=nd_A for all n, 1?n?m/d_A. We also consider the defect index of f(A) for a finite Blaschke product f and show that d_f(A)=d_Aⁿ, where n is the number of zeros of f.     

Determination of a matrix function using the divided difference method of Newton and the interpolation technique of Hermite

Computing a function f(A) of an n-by-n matrix A is a frequently occurring problem in control theory and other applications. In this paper we introduce an effective approach for the determination of matrix function f(A). We propose a new technique which is based on the extension of Newton divided difference and the interpolation technique of Hermite and using the eigenvalues of the given matrix A. The new algorithm is tested on several problems to show the efficiency of the presented method. Finally, the application of this method in control theory is highlighted.     

Google PageRanking problem: The model and the analysis

The spectral and Jordan structures of the Web hyperlink matrix G(c)=cG+(1−c)ev^T have been analyzed when G is the basic (stochastic) Google matrix, c is a real parameter such that 0<c<1, v is a nonnegative probability vector, and e is the all-ones vector. Typical studies have relied heavily on special properties of nonnegative, positive, and stochastic matrices. There is a unique nonnegative vector y(c) such that y(c)^TG(c)=y(c)^T and y(c)^Te=1. This PageRank vector y(c) can be computed effectively by the power method.We consider a square complex matrix A and nonzero complex vectors x and v such that Ax=λx and v^∗x=1. We use standard matrix analytic tools to determine the eigenvalues, the Jordan blocks, and a distinguished left λ-eigenvector of A(c)=cA+(1−c)λxv^∗ as a function of a complex variable c. If λ is a semisimple eigenvalue of A, there is a uniquely determined projection N such that lim_c→1y(c)=Nv for all v; this limit may fail to exist for some v if λ is not semisimple. As a special case of our results, we obtain a complex analog of PageRank for the Web hyperlink matrix G(c) with a complex parameter c. We study regularity, limits, expansions, and conditioning of y(c) and we propose algorithms (e.g., complex extrapolation, power method on a modified matrix etc.) that may provide an efficient way to compute PageRank also with c close or equal to 1. An interpretation of the limit vector Nv and a related critical discussion on the model, on its adherence to reality, and possible ways for its improvement, represent the contribution of the paper on modeling issues.     

Complete hypersurfaces with constant scalar curvature in spheres

To a given immersion ${i:M^n\to \mathbb S^{n+1}}$ with constant scalar curvature R, we associate the supremum of the squared norm of the second fundamental form sup |A|². We prove the existence of a constant C _n (R) depending on R and n so that R ≥ 1 and sup |A|² = C _n (R) imply that the hypersurface is a H(r)-torus ${\mathbb S^1(\sqrt{1-r^2})\times\mathbb S^{n-1} (r)}$ . For R > (n ? 2)/n we use rotation hypersurfaces to show that for each value C > C _n (R) there is a complete hypersurface in ${\mathbb S^{n+1}}$ with constant scalar curvature R and sup |A|² = C, answering questions raised by Q. M. Cheng.