Table of Contents
Fetching ...

A Normal Form Algorithm for Tensor Rank Decomposition

Simon Telen, Nick Vannieuwenhoven

TL;DR

This work addresses the tensor rank decomposition problem by recasting CPD as solving a low-degree polynomial system derived from the kernel of a tensor flattening. It then applies a multi-homogeneous eigenvalue approach via homogeneous normal forms to extract the CPD factors in a non-iterative, linear-algebraic fashion. The authors introduce the notion of multigraded regularity to bound the problem’s complexity and provide conjectures and proven cases that support polynomial-time behavior in many practical tensor formats. Extensive numerical experiments show that the proposed cpd_hnf algorithm achieves higher accuracy and lower memory usage than state-of-the-art methods, even for large, high-order tensors, and demonstrates robustness to noise. The method offers a scalable, algebraically grounded alternative to optimization-based tensor decompositions with potential extensions to higher-order and symmetric tensor settings.

Abstract

We propose a new numerical algorithm for computing the tensor rank decomposition or canonical polyadic decomposition of higher-order tensors subject to a rank and genericity constraint. Reformulating this computational problem as a system of polynomial equations allows us to leverage recent numerical linear algebra tools from computational algebraic geometry. We characterize the complexity of our algorithm in terms of an algebraic property of this polynomial system -- the multigraded regularity. We prove effective bounds for many tensor formats and ranks, which are of independent interest for overconstrained polynomial system solving. Moreover, we conjecture a general formula for the multigraded regularity, yielding a (parameterized) polynomial time complexity for the tensor rank decomposition problem in the considered setting. Our numerical experiments show that our algorithm can outperform state-of-the-art numerical algorithms by an order of magnitude in terms of accuracy, computation time, and memory consumption.

A Normal Form Algorithm for Tensor Rank Decomposition

TL;DR

This work addresses the tensor rank decomposition problem by recasting CPD as solving a low-degree polynomial system derived from the kernel of a tensor flattening. It then applies a multi-homogeneous eigenvalue approach via homogeneous normal forms to extract the CPD factors in a non-iterative, linear-algebraic fashion. The authors introduce the notion of multigraded regularity to bound the problem’s complexity and provide conjectures and proven cases that support polynomial-time behavior in many practical tensor formats. Extensive numerical experiments show that the proposed cpd_hnf algorithm achieves higher accuracy and lower memory usage than state-of-the-art methods, even for large, high-order tensors, and demonstrates robustness to noise. The method offers a scalable, algebraically grounded alternative to optimization-based tensor decompositions with potential extensions to higher-order and symmetric tensor settings.

Abstract

We propose a new numerical algorithm for computing the tensor rank decomposition or canonical polyadic decomposition of higher-order tensors subject to a rank and genericity constraint. Reformulating this computational problem as a system of polynomial equations allows us to leverage recent numerical linear algebra tools from computational algebraic geometry. We characterize the complexity of our algorithm in terms of an algebraic property of this polynomial system -- the multigraded regularity. We prove effective bounds for many tensor formats and ranks, which are of independent interest for overconstrained polynomial system solving. Moreover, we conjecture a general formula for the multigraded regularity, yielding a (parameterized) polynomial time complexity for the tensor rank decomposition problem in the considered setting. Our numerical experiments show that our algorithm can outperform state-of-the-art numerical algorithms by an order of magnitude in terms of accuracy, computation time, and memory consumption.

Paper Structure

This paper contains 22 sections, 19 theorems, 80 equations, 5 figures, 1 table, 2 algorithms.

Key Result

Theorem 1.1

Consider the tensor space $\mathbb{C}^{\ell+1}\otimes\mathbb{C}^{m+1}\otimes\mathbb{C}^{n+1}$ of dimension $M = (\ell+1)(m+1)(n+1)$ with $\ell\ge m \ge n$. If conj:reg holds and $\mathcal{A}$ is a genericA property on a variety $\mathcal{V}$ is "generic" if the locus where the property does not hold

Figures (5)

  • Figure 1: A comparison of the $\log_{10}$ of the relative backward error of the proposed cpd_hnf in the four variants discussed in \ref{['sec_sub_accuracy']} on random rank-$r$ tensors in $\mathbb{R}^{\ell+1} \otimes \mathbb{R}^{m+1} \otimes \mathbb{R}^{n+1}$ with $\ell\ge m\ge n$. The largest dimension and the rank satisfy $\ell+1 = r = \min\{\mathcal{R}(m,n,(2,1)), mn\}$, where the function $\mathcal{R}$ is as in \ref{['eqn_Rc_expression']}. The color scale is the same in all plots.
  • Figure 2: A comparison of the $\log_{10}$ of the relative backward error of the proposed cpd_hnf and the state-of-the-art cpd_ddl from domanov2014canonicaldomanov2017canonical on random rank-$r$ tensors in $\mathbb{R}^{\ell+1} \otimes \mathbb{R}^{m+1} \otimes \mathbb{R}^{n+1}$ with $\ell\ge m\ge n$. The largest dimension and the rank satisfy $\ell+1 = r = \min\{\mathcal{R}(m,n,(d,1)), mn\}$, where the function $\mathcal{R}$ is as in \ref{['eqn_Rc_expression']}. The outcomes for $d = 2, 3, 4$ are shown respectively in the left, middle, and right plots. The color scale is the same in all plots.
  • Figure 3: The $\log_{10}$ of the total execution time (seconds) of cpd_hnf in the setup from \ref{['fig_accuracy']}.
  • Figure 4: The $\log_{10}$ of the size factor $\mu$ of cpd_ddl relative to cpd_hnf. The value $a$ means that cpd_ddl consumes $10^a\times$ the memory cpd_hnf needs to store $R_I(d,e)$.
  • Figure 5: The $\log_{10}$ of the relative backward error of decomposing a random rank-$r$ tensor of size $150 \times 25 \times 10$, corrupted by additive white Gaussian noise of relative magnitude $10^{e}$. The setup is described in detail in \ref{['sec_sub_noisy']}.

Theorems & Definitions (50)

  • Example 1.1: Motivating example
  • Theorem 1.1
  • Lemma 2.1
  • proof
  • Example 2.1
  • Lemma 2.2
  • proof
  • Remark 2.1
  • Example 3.1
  • Lemma 3.1
  • ...and 40 more