On the perturbation of rank‐one symmetric tensors

The problem of symmetric rank‐one approximation of symmetric tensors is important in independent components analysis, also known as blind source separation, as well as polynomial optimization. We derive several perturbative results that are relevant to the well‐posedness of recovering rank‐one structure from approximately‐rank‐one symmetric tensors. We also specialize the analysis of the shifted symmetric higher‐order power method, an algorithm for computing symmetric tensor eigenvectors, to approximately‐rank‐one symmetric tensors. Copyright © 2012 John Wiley & Sons, Ltd.     

A preconditioned low‐rank CG method for parameter‐dependent Lyapunov matrix equations

This paper is concerned with the numerical solution of symmetric large‐scale Lyapunov equations with low‐rank right‐hand sides and coefficient matrices depending on a parameter. Specifically, we consider the situation when the parameter dependence is sufficiently smooth, and the aim is to compute solutions for many different parameter samples. On the basis of existing results for Lyapunov equations and parameter‐dependent linear systems, we prove that the tensor containing all solution samples typically allows for an excellent low multilinear rank approximation. Stacking all sampled equations into one huge linear system, this fact can be exploited by combining the preconditioned CG method with low‐rank truncation. Our approach is flexible enough to allow for a variety of preconditioners based, for example, on the sign function iteration or the alternating direction implicit method. Copyright © 2013 John Wiley & Sons, Ltd.     

Low‐rank approximation of tensors via sparse optimization

The goal of this paper is to find a low‐rank approximation for a given nth tensor. Specifically, we give a computable strategy on calculating the rank of a given tensor, based on approximating the solution to an NP‐hard problem. In this paper, we formulate a sparse optimization problem via an l₁‐regularization to find a low‐rank approximation of tensors. To solve this sparse optimization problem, we propose a rescaling algorithm of the proximal alternating minimization and study the theoretical convergence of this algorithm. Furthermore, we discuss the probabilistic consistency of the sparsity result and suggest a way to choose the regularization parameter for practical computation. In the simulation experiments, the performance of our algorithm supports that our method provides an efficient estimate on the number of rank‐one tensor components in a given tensor. Moreover, this algorithm is also applied to surveillance videos for low‐rank approximation.     

Quantized CP approximation and sparse tensor interpolation of function‐generated data

In this article, we consider the iterative schemes to compute the canonical polyadic (CP) approximation of quantized data generated by a function discretized on a large uniform grid in an interval on the real line. This paper continues the research on the quantics‐tensor train (QTT) method (“O(d log N)‐quantics approximation of N‐d tensors in high‐dimensional numerical modeling” in Constructive Approximation, 2011) developed for the tensor train (TT) approximation of the quantized images of function related data. In the QTT approach, the target vector of length 2^L is reshaped to a Lth‐order tensor with two entries in each mode (quantized representation) and then approximated by the QTT tensor including 2r²L parameters, where r is the maximal TT rank. In what follows, we consider the alternating least squares (ALS) iterative scheme to compute the rank‐r CP approximation of the quantized vectors, which requires only 2rL?2^L parameters for storage. In the earlier papers (“Tensors‐structured numerical methods in scientific computing: survey on recent advances” in Chemom Intell Lab Syst, 2012), such a representation was called Q_Can format, whereas in this paper, we abbreviate it as the QCP (quantized canonical polyadic) representation. We test the ALS algorithm to calculate the QCP approximation on various functions, and in all cases, we observed the exponential error decay in the QCP rank. The main idea for recovering a discretized function in the rank‐r QCP format using the reduced number of the functional samples, calculated only at O(2rL) grid points, is presented. The special version of the ALS scheme for solving the arising minimization problem is described. This approach can be viewed as the sparse QCP‐interpolation method that allows to recover all 2rL representation parameters of the rank‐r QCP tensor. Numerical examples show the efficiency of the QCP‐ALS‐type iteration and indicate the exponential convergence rate in r.     

A sequential subspace projection method for extreme Z‐eigenvalues of supersymmetric tensors

Z‐eigenvalues of tensors, especially extreme ones, are quite useful and are related to many problems, such as automatic control, quantum physics, and independent component analysis. For supersymmetric tensors, calculating the smallest/largest Z‐eigenvalue is equivalent to solving a global minimization/maximization problem of a homogenous polynomial over the unit sphere. In this paper, we utilize the sequential subspace projection method (SSPM) to find extreme Z‐eigenvalues and the corresponding Z‐eigenvectors. The main idea of SSPM is to form a 2‐dimensional subspace at the current point and then solve the original optimization problem in the subspace. SSPM benefits from the fact that the 2‐dimensional subproblem can be solved by a direct method. Global convergence and linear convergence are established for supersymmetric tensors under certain assumptions. Preliminary numerical results over several testing problems show that SSPM is very promising. Besides, the globalization strategy of random phase can be easily incorporated into SSPM, which promotes the ability to find extreme Z‐eigenvalues. Copyright © 2014 John Wiley & Sons, Ltd.     

A constructive arbitrary‐degree Kronecker product decomposition of tensors

We generalize the matrix Kronecker product to tensors and propose the tensor Kronecker product singular value decomposition that decomposes a real k‐way tensor into a linear combination of tensor Kronecker products with an arbitrary number of d factors. We show how to construct , where each factor is also a k‐way tensor, thus including matrices (k=2) as a special case. This problem is readily solved by reshaping and permuting into a d‐way tensor, followed by a orthogonal polyadic decomposition. Moreover, we introduce the new notion of general symmetric tensors (encompassing symmetric, persymmetric, centrosymmetric, Toeplitz and Hankel tensors, etc.) and prove that when is structured then its factors will also inherit this structure.     

The tensor t‐function: A definition for functions of third‐order tensors

A definition for functions of multidimensional arrays is presented. The definition is valid for third‐order tensors in the tensor t‐product formalism, which regards third‐order tensors as block circulant matrices. The tensor function definition is shown to have similar properties as standard matrix function definitions in fundamental scenarios. To demonstrate the definition's potential in applications, the notion of network communicability is generalized to third‐order tensors and computed for a small‐scale example via block Krylov subspace methods for matrix functions. A complexity analysis for these methods in the context of tensors is also provided.     

Decompositions of third‐order tensors: HOSVD,T‐SVD,and Beyond

The higher order singular value decomposition, which is regarded as a generalization of the matrix singular value decomposition (SVD), has a long history and is well established, while the T‐SVD is relatively new and lacks systematic analysis. Because of the unusual tensor‐tensor product that the T‐SVD is based on, the form of the T‐SVD may be difficult to comprehend. The main aim of this article is to establish a connection between these two decompositions. By converting the form of the T‐SVD into the sum of outer product terms, we compare the forms of the two decompositions. Moreover, from establishing the connection, a new decomposition which has a specific nonzero pattern, is proposed and developed. Numerical examples are given to demonstrate the useful ability of the new decomposition for approximation and data compression.     

Fast tensor method for summation of long‐range potentials on 3D lattices with defects

In this paper, we present a method for fast summation of long‐range potentials on 3D lattices with multiple defects and having non‐rectangular geometries, based on rank‐structured tensor representations. This is a significant generalization of our recent technique for the grid‐based summation of electrostatic potentials on the rectangular L × L × L lattices by using the canonical tensor decompositions and yielding the O(L) computational complexity instead of O(L³) by traditional approaches. The resulting lattice sum is calculated as a Tucker or canonical representation whose directional vectors are assembled by the 1D summation of the generating vectors for the shifted reference tensor, once precomputed on large N × N × N representation grid in a 3D bounding box. The tensor numerical treatment of defects is performed in an algebraic way by simple summation of tensors in the canonical or Tucker formats. To diminish the considerable increase in the tensor rank of the resulting potential sum, the ?‐rank reduction procedure is applied based on the generalized reduced higher‐order SVD scheme. For the reduced higher‐order SVD approximation to a sum of canonical/Tucker tensors, we prove the stable error bounds in the relative norm in terms of discarded singular values of the side matrices. The required storage scales linearly in the 1D grid‐size, O(N), while the numerical cost is estimated by O(NL). The approach applies to a general class of kernel functions including those for the Newton, Slater, Yukawa, Lennard‐Jones, and dipole‐dipole interactions. Numerical tests confirm the efficiency of the presented tensor summation method; we demonstrate that a sum of millions of Newton kernels on a 3D lattice with defects/impurities can be computed in seconds in Matlab implementation. The tensor approach is advantageous in further functional calculus with the lattice potential sums represented on a 3D grid, like integration or differentiation, using tensor arithmetics of 1D complexity. Copyright © 2015 John Wiley & Sons, Ltd.     

CP decomposition for tensors via alternating least squares with QR decomposition

The CP tensor decomposition is used in applications such as machine learning and signal processing to discover latent low-rank structure in multidimensional data. Computing a CP decomposition via an alternating least squares (ALS) method reduces the problem to several linear least squares problems. The standard way to solve these linear least squares subproblems is to use the normal equations, which inherit special tensor structure that can be exploited for computational efficiency. However, the normal equations are sensitive to numerical ill-conditioning, which can compromise the results of the decomposition. In this paper, we develop versions of the CP-ALS algorithm using the QR decomposition and the singular value decomposition, which are more numerically stable than the normal equations, to solve the linear least squares problems. Our algorithms utilize the tensor structure of the CP-ALS subproblems efficiently, have the same complexity as the standard CP-ALS algorithm when the input is dense and the rank is small, and are shown via examples to produce more stable results when ill-conditioning is present. Our MATLAB implementation achieves the same running time as the standard algorithm for small ranks, and we show that the new methods can obtain lower approximation error.     

A semidefinite program approach for computing the maximum eigenvalue of a class of structured tensors and its applications in hypergraphs and copositivity test

Finding the maximum eigenvalue of a symmetric tensor is an important topic in tensor computation and numerical multilinear algebra. In this paper, we introduce a new class of structured tensors called W‐tensors, which not only extends the well‐studied nonnegative tensors by allowing negative entries but also covers several important tensors arising naturally from spectral hypergraph theory. We then show that finding the maximum H‐eigenvalue of an even‐order symmetric W‐tensor is equivalent to solving a structured semidefinite program and hence can be validated in polynomial time. This yields a highly efficient semidefinite program algorithm for computing the maximum H‐eigenvalue of W‐tensors and is based on a new structured sums‐of‐squares decomposition result for a nonnegative polynomial induced by W‐tensors. Numerical experiments illustrate that the proposed algorithm can successfully find the maximum H‐eigenvalue of W‐tensors with dimension up to 10,000, subject to machine precision. As applications, we provide a polynomial time algorithm for computing the maximum H‐eigenvalues of large‐size Laplacian tensors of hyperstars and hypertrees, where the algorithm can be up to 13 times faster than the state‐of‐the‐art numerical method introduced by Ng, Qi, and Zhou in 2009. Finally, we also show that the proposed algorithm can be used to test the copositivity of a multivariate form associated with symmetric extended Z‐tensors, whose order may be even or odd.     

A nonlinearly preconditioned conjugate gradient algorithm for rank‐R canonical tensor approximation

Alternating least squares (ALS) is often considered the workhorse algorithm for computing the rank‐R canonical tensor approximation, but for certain problems, its convergence can be very slow. The nonlinear conjugate gradient (NCG) method was recently proposed as an alternative to ALS, but the results indicated that NCG is usually not faster than ALS. To improve the convergence speed of NCG, we consider a nonlinearly preconditioned NCG (PNCG) algorithm for computing the rank‐R canonical tensor decomposition. Our approach uses ALS as a nonlinear preconditioner in the NCG algorithm. Alternatively, NCG can be viewed as an acceleration process for ALS. We demonstrate numerically that the convergence acceleration mechanism in PNCG often leads to important pay‐offs for difficult tensor decomposition problems, with convergence that is significantly faster and more robust than for the stand‐alone NCG or ALS algorithms. We consider several approaches for incorporating the nonlinear preconditioner into the NCG algorithm that have been described in the literature previously and have met with success in certain application areas. However, it appears that the nonlinearly PNCG approach has received relatively little attention in the broader community and remains underexplored both theoretically and experimentally. Thus, this paper serves several additional functions, by providing in one place a concise overview of several PNCG variants and their properties that have only been described in a few places scattered throughout the literature, by systematically comparing the performance of these PNCG variants for the tensor decomposition problem, and by drawing further attention to the usefulness of nonlinearly PNCG as a general tool. In addition, we briefly discuss the convergence of the PNCG algorithm. In particular, we obtain a new convergence result for one of the PNCG variants under suitable conditions, building on known convergence results for non‐preconditioned NCG. Copyright © 2014 John Wiley & Sons, Ltd.     

Shifted power method for computing tensor H‐eigenpairs

In this paper, we propose a shifted symmetric higher‐order power method for computing the H‐eigenpairs of a real symmetric even‐order tensor. The local convergence of the method is proved. In addition, by utilizing the fixed‐point analysis, we can characterize exactly which H‐eigenpairs can be found and which cannot be found by the method. Numerical examples are presented to illustrate the performance of the method. Copyright © 2016 John Wiley & Sons, Ltd.     

News Algorithms for tensor decomposition based on a reduced functional

We study the least squares functional of the canonical polyadic tensor decomposition for third order tensors by eliminating one factor matrix, which leads to a reduced functional. An analysis of the reduced functional leads to several equivalent optimization problem, such as a Rayleigh quotient or a projection. These formulations are the basis of several new algorithms as follows: the Centroid Projection method for efficient computation of suboptimal solutions and fixed‐point iteration methods for approximating the best rank‐1 and the best rank‐R decompositions under certain nondegeneracy conditions. Copyright © 2013 John Wiley & Sons, Ltd.     

Reeb parallel Ricci tensor for homogeneous real hypersurfaces in complex hyperbolic two‐plane Grassmannians

In this paper, we introduce the notion of Reeb parallel Ricci tensor for homogeneous real hypersurfaces in complex hyperbolic two‐plane Grassmannians which has a remarkable geometric structure as a Hermitian symmetric space of rank 2. By using a new method of simultaneous diagonalizations, we give a complete classification for real hypersurfaces in complex hyperbolic two‐plane Grassmannians with the Reeb parallel Ricci tensor.     

A new two‐phase structure‐preserving doubling algorithm for critically singular M‐matrix algebraic Riccati equations

Among numerous iterative methods for solving the minimal nonnegative solution of an M‐matrix algebraic Riccati equation, the structure‐preserving doubling algorithm (SDA) stands out owing to its overall efficiency as well as accuracy. SDA is globally convergent and its convergence is quadratic, except for the critical case for which it converges linearly with the linear rate 1/2. In this paper, we first undertake a delineatory convergence analysis that reveals that the approximations by SDA can be decomposed into two components: the stable component that converges quadratically and the rank‐one component that converges linearly with the linear rate 1/2. Our analysis also shows that as soon as the stable component is fully converged, the rank‐one component can be accurately recovered. We then propose an efficient hybrid method, called the two‐phase SDA, for which the SDA iteration is stopped as soon as it is determined that the stable component is fully converged. Therefore, this two‐phase SDA saves those SDA iterative steps that previously have to have for the rank‐one component to be computed accurately, and thus essentially, it can be regarded as a quadratically convergent method. Numerical results confirm our analysis and demonstrate the efficiency of the new two‐phase SDA. Copyright © 2015 John Wiley & Sons, Ltd.     

On algorithms for and computing with the tensor ring decomposition

Tensor decompositions such as the canonical format and the tensor train format have been widely utilized to reduce storage costs and operational complexities for high‐dimensional data, achieving linear scaling with the input dimension instead of exponential scaling. In this paper, we investigate even lower storage‐cost representations in the tensor ring format, which is an extension of the tensor train format with variable end‐ranks. Firstly, we introduce two algorithms for converting a tensor in full format to tensor ring format with low storage cost. Secondly, we detail a rounding operation for tensor rings and show how this requires new definitions of common linear algebra operations in the format to obtain storage‐cost savings. Lastly, we introduce algorithms for transforming the graph structure of graph‐based tensor formats, with orders of magnitude lower complexity than existing literature. The efficiency of all algorithms is demonstrated on a number of numerical examples, and in certain cases, we demonstrate significantly higher compression ratios when compared to previous approaches to using the tensor ring format.     

Convergence analysis of projected fixed‐point iteration on a low‐rank matrix manifold

In this paper, we analyze the convergence of a projected fixed‐point iteration on a Riemannian manifold of matrices with fixed rank. As a retraction method, we use the projector splitting scheme. We prove that the convergence rate of the projector splitting scheme is bounded by the convergence rate of standard fixed‐point iteration without rank constraints multiplied by the function of initial approximation. We also provide counterexample to the case when conditions of the theorem do not hold. Finally, we support our theoretical results with numerical experiments.     

A practical method for computing the largest M‐eigenvalue of a fourth‐order partially symmetric tensor

In this paper, we consider a bi‐quadratic homogeneous polynomial optimization problem over two unit spheres arising in nonlinear elastic material analysis and in entanglement studies in quantum physics. The problem is equivalent to computing the largest M‐eigenvalue of a fourth‐order tensor. To solve the problem, we propose a practical method whose validity is guaranteed theoretically. To make the sequence generated by the method converge to a good solution of the problem, we also develop an initialization scheme. The given numerical experiments show the effectiveness of the proposed method. Copyright © 2009 John Wiley & Sons, Ltd.     

Convergence of a transition probability tensor of a higher‐order Markov chain to the stationary probability vector

In this paper, first we introduce a new tensor product for a transition probability tensor originating from a higher‐order Markov chain. Subsequently, some properties of the new tensor product are explained, and its relationship with the stationary probability vector is studied. Also, similarity between results obtained by this new product and the first‐order case is shown. Furthermore, we prove the convergence of a transition probability tensor to the stationary probability vector. Finally, we show how to achieve a stationary probability vector with some numerical examples and make some comparison between the proposed method and another existing method for obtaining stationary probability vectors. Copyright © 2016 John Wiley & Sons, Ltd.