# Generalized $\mathcal{L}$-product for hight order tensors with applications
using GPU computations.

A.H Bentbib Faculté des Sciences et Techniques-Gueliz, Laboratoire de
Mathématiques Appliquées et Informatique, Marrakech, Morocco M.
Elalj44footnotemark: 4 A. El Hachimi44footnotemark: 4 K. Jbilou44footnotemark:
4 Université du Littoral Cote d’Opale, LMPA, 50 rue F. Buisson, 62228 Calais-
Cedex, France  A. Ratnani laboratory MSDA, Mohammed VI Polytechnic University,
Green City, Morocco

###### Abstract

In this paper, we will present a generalization of the $\mathcal{L}$-tensor
product ($*_{\mathcal{L}}$-product) including generalization of the well known
tensor cosine and T-products that were defined for third-order tensors and
based on fast Fourier transform and discrete cosine transform (DCT). We will
give some applications on tensor completion. To solve some optimization
problems linked with the problem of tensor completion, we will use the
Proximal Gradient Algorithm (PGA) to solve some derived optimization problems.
Numerical tests are given to show the effectiveness of the proposed methods
and also present some tests using GPU computation.

###### keywords:

Tensor completion, tensor nuclear norm, tensor $*_{\mathcal{L}}$-product,
tensor $*_{\mathcal{L}}$-svd, GPU.

##  1 Introduction

In the last decade, tensors become an important multilinear algebra tool
involved in many modern problems such completion [9, 18, 24], principal
component analysis [13], image processing [20, 12, 4] and others.
The classical $n$-mode product leads to many concepts and developements when
working with multidimensional data. The CP and the Tucker compressions were
introduced as natural generalization of the classical singular value
decomposition (SVD) for matrices; see [15, 4, 12, 13, 24].
In the last years, new tensor-tensor products such as cosine-product
(c-product), using discrete cosine or T-product, using Fast Fourier Transform
(FFT), were introduced for third-order tensors, studied and applied to image
processing and other fields; see [19, 1, 25, 13, 20]. In the present paper, we
generalize those tensor-tensor products for high-order tensors. Using those
new products, we will propose new completion models. We give some theoretical
results and some numerical examples in color video processing.
The outline of this paper will as follows In Section 2 we will give some
definitions and remind some known results of the third-order tensor-tensor
product based on Fast Fourier Transform and cosine transform. In Section 3, we
present our new generalized c-product for any order of tensors and give some
important results. Section 4 presents a novel models of tensor completion for
tensors of any order by using the PGA. Section 5 will be devoted to some
numerical experiments with someexperiments using GPU.

##  2 Definitions and notations

In this subsection we will define some notions that will help us in the paper.
We denote tensors by Euler script letters, e.g., $\mathcal{X}$, matrices will
be denoted by boldface capital letters ,e.g., $X$, vectors by boldface
lowercase letters e.g., x and scalars by lowercase letters, e.g., $x$. Also we
will denote the $(i_{1},\,i_{2},\,...,\,i_{K})^{th}$ for a $K^{th}$-order
tensor $\mathcal{X}$ by $\mathcal{X}_{i_{1},\,i_{2},\,...,\,i_{K}}$. Also we
will denote $\mathbb{C}^{I_{1}\times I_{2}\times\dots\times I_{K}}$ by
$\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}\times\dots\times I_{K}}$, and the
$K^{th}$-order tensor-scalar space as $\mathbb{K}_{I_{3}\times\dots\times
I_{K}}=\mathbb{C}^{1\times 1\times I_{3}\times\dots\times I_{K}}$.
Let $\mathcal{A}$ and $\mathcal{B}$ be two $K^{th}$-order tensors in
$\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}\times\dots\times I_{k}}$, we will
define the inner product between $\mathcal{A}$ and $\mathcal{B}$ by

$\left<\mathcal{A},\,\mathcal{B}\right>=\sum_{i_{1},\,i_{2},\,\dots,\,i_{K}=1}^{I_{1},\,I_{2},\,\dots,\,I_{K}}\mathcal{A}_{i_{1},\,i_{2},\,\dots,\,i_{K}}\mathcal{B}_{i_{1},\,i_{2},\,\dots,\,i_{K}}.$
(1)

The associated norm is defined by

$\left\|\mathcal{A}\right\|_{F}^{2}=\sum_{i_{1},\,i_{2},\,\dots,\,i_{K}=1}^{I_{1},\,I_{2},\,\dots,\,I_{K}}\mathcal{A}_{i_{1},\,i_{2},\,\dots,\,i_{K}}^{2}.$
(2)

The notion of columns and rows of matrices are extended to the case of
tensors, where in this case we talk about $n$-mode fiber instead of columns
and rows, with the $n$-mode fiber of a $K^{th}$-order tensor $\mathcal{X}$ is
defined by fixing all the indexes except the $n^{th}$ one.
There is some ways to transform a tensor to a matrix which consist to make the
operations on tensors easier, there is for example the $n$-mode matricization
[4, 24] defined as follows

######  Definition 1.

Let $\mathcal{X}\in\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}\times\dots\times
I_{N}}$, then the $n$-mode matricization of $\mathcal{X}$ denoted by
$X_{(n)}\in\mathbb{K}^{I_{n}\times\prod_{k=1,\,k\neq n}^{N}I_{k}}$ and it is
defined by making the $n$-mode fibers as columns of $X_{(n)}$, i.e., the
$(i_{1},\,i_{2},\,\dots,\,i_{N})^{th}$ element of $\mathcal{X}$ maps to a
matrix element $(i_{n},\,j)$ satisfying

$j=1+\sum_{k=1,\,j\neq n}^{N}(i_{k}-1)J_{k}\,\,for\,\,J_{k}=\prod_{m=1,\,m\neq
n}^{N}I_{m}.$

For third-order tensors $\mathcal{X}\in\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}}$, $X_{(1)}$, $X_{(2)}$ and $X_{(3)}$ are given by

$\displaystyle X_{(1)}$ $\displaystyle=$
$\displaystyle\left[\mathcal{X}^{(1)},\,\mathcal{X}^{(2)},\,\dots,\,\mathcal{X}^{(I_{3})}\right],$
$\displaystyle X_{(2)}$ $\displaystyle=$
$\displaystyle\left[\left(\mathcal{X}^{(1)}\right)^{T},\,\left(\mathcal{X}^{(2)}\right)^{T},\,\dots,\,\left(\mathcal{X}^{(I_{3})}\right)^{T}\right],$
$\displaystyle X_{(3)}$ $\displaystyle=$ $\displaystyle\left[{\tt
sq}\left(\mathcal{X}(:,\,1,\,:)\right)^{T},\,{\tt
sq}\left(\mathcal{X}(:,\,2,\,:)\right)^{T},\,\dots,\,{\tt
sq}\left(\mathcal{X}(:,\,I_{2},\,:)\right)^{T}\right],$

where $\mathcal{X}^{(i)}$ denotes the $i^{th}$ frontal slice
$\left(\mathcal{X}^{(i)}=\mathcal{X}(:,:,i)\right)$ and ${\tt sq}$ transforms
the tensor $\mathcal{X}\in\mathbb{K}^{I_{1}\times 1}_{I_{3}}$ to a matrix
$X\in\mathbb{K}^{I_{1}\times I_{3}}$, i.e., $X={\tt
sq}\left(\mathcal{X}\right)$.
The $n$-mode product, which is a product between a tensor and a matrix in the
$n$-mode [4] is defined in the following definition

######  Definition 2.

Let $\mathcal{X}\in\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}\times\dots\times
I_{K}}$ and $U\in\mathbb{K}^{J\times I_{n}}$ where $J$ is a positive nonzero
integer. Then the $n$-mode product $\mathcal{X}\times_{n}U$ is the tensor in
$\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}\times\dots\times I_{n-1}\times J\times
I_{n+1}\times\dots\times I_{K}}$, where its
$\left(i_{1},\dots,i_{n-1},j,i_{n+1},\dots,i_{K}\right)^{th}$ element is
defined by

$\left(\mathcal{X}\times_{n}U\right)_{i_{1},\,\dots,\,i_{n-1},\,j,\,i_{n+1},\,\dots,\,i_{K}}=\sum_{i_{n}=1}^{I_{n}}\mathcal{X}_{i_{1},\,i_{2},\,\dots,\,i_{K}}U_{j,\,i_{n}}.$
(3)

Some useful properties of the $n$-mode product are given as follows. Let the
tensor $\mathcal{X}$, and the matrices $U$ and $V$ of appropriate sizes, then

$\displaystyle\mathcal{Y}=\mathcal{X}\times_{n}U$
$\displaystyle\Longleftrightarrow$ $\displaystyle Y_{(n)}=UX_{(n)},\;and$
$\displaystyle\mathcal{X}\times_{n}U\times_{m}V$ $\displaystyle=$
$\displaystyle\mathcal{X}\times_{m}V\times_{n}U.$

We will also use the notion of tensor face-wise product defined next.

######  Definition 3.

Let $\mathcal{X}\in\mathbb{K}^{I_{1}\times n}_{I_{3}}$ and
$\mathcal{Y}\in\mathbb{K}^{n\times I_{2}}_{I_{3}}$ two third-order tensors,
then the face-wise product between $\mathcal{X}$ and $\mathcal{Y}$ is given by
the tensor of size $I_{1}\times I_{2}\times I_{3}$ where its $i^{th}$ frontal
slice is given from the product between the $i^{th}$ frontal slices of
$\mathcal{X}$ and $\mathcal{Y}$, i.e.,

$\left(\mathcal{X}\triangle\mathcal{Y}\right)^{(i)}=\mathcal{X}^{(i)}\triangle\mathcal{Y}^{(i)}.$
(4)

Classical tensor decompositions such as CP decomposition [4], Tucker
decomposition [24], block term decomposition [16] give nice results in many
tensor applications. However, those decompositions suffer from the high
computational cost for large problems. In the recent years new tensor
decompositions of the third-order case and based on tensor-tensor product
using the Fourier domain such as the t-product [12] and cosine-product
(c-product) [2], have been defined and used for many image processing
applications; see [1, 3, 9, 13].
In this section we will try to remind the most important results of those
types of tensor-tensor product.
The main idea of this type of tensor products is to transform the tensors to
another domain which is called the transform domain, like Fourier domain,
cosine domain. Then all the operations are done in the transformed domain
using for examlple FFT on each tube to speed-up the executing time. This kind
of transformation could be defined in the following.

######  Definition 4.

Let $M$ be an invertible matrix of size $I_{3}\times I_{3}$, we define the
operator $\mathcal{L}$ as

$\displaystyle\mathcal{L}:\;\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}}$
$\displaystyle\longrightarrow$ $\displaystyle\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}}$ $\displaystyle\mathcal{A}$ $\displaystyle\longrightarrow$
$\displaystyle\mathcal{A}\times_{3}M$

and its inverse is defined as

$\displaystyle\mathcal{L}^{-1}:\;\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}}$
$\displaystyle\longrightarrow$ $\displaystyle\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}}$ $\displaystyle\mathcal{A}$ $\displaystyle\longrightarrow$
$\displaystyle\mathcal{A}\times_{3}M^{-1}$

Now we can define the tensor-tensor product of two third-order tensors.

######  Definition 5.

Let $\mathcal{L}$ be an invertible operator, then the tensor-tensor product
between two third-order tensors associated with the operator $\mathcal{L}$, is
denoted by $*_{\mathcal{L}}$ and is given by

$\mathcal{A}*_{\mathcal{L}}\mathcal{B}=\mathcal{L}^{-1}\left(\mathcal{L}\left(\mathcal{A}\right)\triangle\mathcal{L}\left(\mathcal{B}\right)\right).$
(5)

where the tensors $\mathcal{A}$ in $\mathbb{K}^{I_{1}\times l}_{I_{3}}$ and
$\mathcal{B}$ in $\mathbb{K}^{l\times I_{2}}_{I_{3}}$.

The matrix $M$ in Definition 4 depends on the type of the product, for example
if we use the t-product [12], the matrix $M$ is the matrix of discrete Fourier
transform $F_{I_{3}}$ where the Fourier matrix $F_{n}\in\mathbb{C}^{n\times
n}$ is given by

$F_{n}=[\omega_{n}^{(i-1)(j-1)}];\,i,j=1,\ldots,n-1;\;and\;\omega_{n}=e^{\frac{-2i\pi}{n}}.$
(6)

for $n\in\mathbb{N}^{*}$. Notice that $\dfrac{F_{n}}{\sqrt{n}}$ is unitary,
i.e., $F_{n}F_{n}^{*}=nI_{n}$.
In the case of c-product [2], the matrix $M$ is defined as

$M=W^{-1}_{I_{3}}C_{I_{3}}\left(I_{I_{3}}+Z_{I_{3}}\right),$ (7)

where $W_{I_{3}}={\tt diag}\left(C_{I_{3}}(:,1)\right)$, the matrix
$Z_{I_{3}}$ is the circulant upshift matrix defined by $Z_{I_{3}}={\tt
diag(ones(I_{3}-1,1),1)}$ and $C_{I_{3}}$ is the matrix of discrete cosine
transform of size $I_{3}\times I_{3}$ and its $(i,j)^{th}$ element is defined
as

$\left(C_{I_{3}}\right)_{i,\,j}=\sqrt{\dfrac{2-\delta_{i,\,j}}{I_{3}}}\cos\left(\dfrac{(i-1)(2j-1)\pi}{2I_{3}}\right)\;\;with\;\;1\leq
i,\,j\leq I_{3},$ (8)

where $\delta_{i,\,j}$ is the Kronecker symbol. Notice that the matrix $C_{n}$
is orthogonal for all $n\in\mathbb{N}^{*}$. We have also to mention that, in
this case also, the matrix $M$ is invertible and
$M^{-1}_{I_{3}}=\left(I_{I_{3}}+Z_{I_{3}}\right)^{-1}C^{*}_{I_{3}}W_{I_{3}}.$
Using those tensor-tensor products, all the classical matrix decomposition,
such as svd, QR and Shur decompositions have been generalized to the tensor
case; see [2, 7]. Many applications of tensor-tensor product use some
optimization algorithms and in our present work, we will use the Proximal
Gradient Algorithm [10] in tensor completion. The method consists in solving
the optimization problem

$\displaystyle\underset{\mathcal{X}\in\mathbb{H}}{\min}$
$\displaystyle\,g(\mathcal{X})$ (9) $\displaystyle{\tt s.t}$
$\displaystyle\textbf{A}\left(\mathcal{X}\right)=\mathcal{B},$

where $\mathbb{H}$ is an Hilbert space equipped with a norm
$\left\|.\right\|$, $g$ is a continuous function, A is a linear map and
$\mathcal{B}$ is an observation.
By referring to [5, 14], this optimization problem can be solved by solving
the following one

$\underset{\mathcal{X}\in\mathbb{H}}{\min}\,\\{F(\mathcal{X})=\mu
g(\mathcal{X})+f(\mathcal{X})\\},$ (10)

where
$f(\mathcal{X})=\dfrac{1}{2}\left\|\mathcal{A}(\mathcal{X})-\mathcal{B}\right\|^{2}$
and $\mu>0$ is the relaxation parameter. The penality function $f$ is convex
and smooth with Lipshitz continuous gradient, with Lipshitz constant $l_{f}$.
To solve (10), we minimize the quadratic function
$Q(\mathcal{X},\mathcal{Y})$, where $\mathcal{Y}$ is chosen and $Q$ is defined
as follows

$Q\left(\mathcal{X},\mathcal{Y}\right)=\mu
g(\mathcal{X})+f(\mathcal{Y})+\left<\nabla
f(\mathcal{Y}),\,\mathcal{X}-\mathcal{Y}\right>+\dfrac{l_{f}}{2}\left\|\mathcal{X}-\mathcal{Y}\right\|^{2}.$
(11)

Getting a solution $\mathcal{X}$ satisfying (11) is equivalent to solve the
following minimization problem

$\underset{\mathcal{X}\in\mathbb{H}}{\min}\,Q(\mathcal{X},\mathcal{Y})=\underset{\mathcal{X}\in\mathbb{H}}{\min}\,\mu
g(\mathcal{X})+\dfrac{l_{f}}{2}\left\|\mathcal{X}-\mathcal{G}\right\|^{2},$
(12)

where $\mathcal{G}=\mathcal{Y}-\dfrac{1}{l_{f}}\nabla f(\mathcal{Y})$. The
problem is solved iteratively by computing $\mathcal{X}^{p+1}$ such that

$\mathcal{X}^{p+1}=\underset{\mathcal{X}\in\mathbb{H}}{\arg\,\min}\,Q(\mathcal{X},\,\mathcal{Y}^{p}),$
(13)

In [14], $\mathcal{Y}^{p}$ was computed by
$\mathcal{Y}^{p}=\mathcal{X}^{p}+\dfrac{t_{p+1}-1}{t_{p}}\left(\mathcal{X}^{p}-\mathcal{X}^{p-1}\right)$
instead of $\mathcal{Y}^{p}=\mathcal{X}^{p}$ for computationally reasons and
$t_{p+1}=\dfrac{1+\sqrt{4t_{p}^{2}+1}}{2}$. The steeps of this algorithm can
be summarized in the following algorithm

Algorithm 1 Proximal Gradient Algorithm (PGA).

1: while not converged do

2:
$\mathcal{Y}^{p}=\mathcal{X}^{p}+\dfrac{t_{p-1}-1}{t_{p}}\left(\mathcal{X}^{p}-\mathcal{X}^{p-1}\right)$.

3: $\mathcal{G}^{p}=\mathcal{Y}^{p}-\dfrac{1}{l_{f}}\nabla
f\left(\mathcal{Y}^{p}\right)$.

4: $\mathcal{X}^{p+1}=\underset{\mathcal{X}\in\mathbb{H}}{\arg\,\min}\,\mu
g(\mathcal{X})+\dfrac{l_{f}}{2}\left\|\mathcal{X}-\mathcal{G}^{p}\right\|^{2}$.

5: $t_{p+1}=\dfrac{1+\sqrt{4t_{p}^{2}+1}}{2}$.

6: $p=p+1.$

7: end while

##  3 Generalized tensor-tensor cosine product

The main inconvenient of tensor-tensor products above is the fact that they
could be used only for third-order tensors. In [27], the authors proposed a
generalization of the t-product and in the present work we propse the
generalization of the c-product with some applications. We will first recall
some important results linked with the c-product for third-order tensors
described in [2]. For two third-order tensors
$\mathcal{A}\in\mathbb{K}^{I_{1}\times l}_{I_{3}}$ and
$\mathcal{B}\in\mathbb{K}^{l\times I_{2}}_{I_{3}}$, the c-product
$\mathcal{A}*_{c}\mathcal{B}$ is defined by

$\mathcal{A}*_{c}\mathcal{B}={\tt ten}\left({\tt
btph}\left(\mathcal{A}\right){\tt
btph}\left(\mathcal{B}\right)\right)\in\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}},$ (14)

where ${\tt btph}$ represents the block-Toeplitz-plus-Hankel matrix defined as

${\tt
btph}\left(\mathcal{A}\right)=\begin{pmatrix}\mathcal{A}^{(1)}&\dots&\dots&\mathcal{A}^{(I_{3})}\\\
\mathcal{A}^{(2)}&\mathcal{A}^{(1)}&\dots&\mathcal{A}^{(I_{3}-1)}\\\
\vdots&\ddots&\ddots&\vdots\\\
\mathcal{A}^{(I_{3})}&\dots&\dots&\mathcal{A}^{(1)}\end{pmatrix}+\begin{pmatrix}\mathcal{A}^{(2)}&\dots&\mathcal{A}^{(I_{3})}&\textbf{0}\\\
\vdots&\reflectbox{$\ddots$}&\reflectbox{$\ddots$}&\mathcal{A}^{(I_{3})}\\\
\mathcal{A}^{(I_{3})}&\textbf{0}&\reflectbox{$\ddots$}&\vdots\\\
\textbf{0}&\mathcal{A}^{(I_{3})}&\dots&\mathcal{A}^{(2)}\end{pmatrix}\in\mathbb{K}^{I_{1}I_{3}\times
I_{2}l}$ (15)

and the operator ${\tt ten}$ is the operator to get back a tensor from its
${\tt btph}$, where ${\tt ten}\left({\tt
btph}\left(\mathcal{A}\right)\right)=\mathcal{A}.$ Computing this product by
the above formula can be expensive, because the matrix ${\tt btph}$ may be
very large. Some properties of ${\tt btph}$ are given in [2]. Among them, the
fact that for a tensor $\mathcal{A}\in\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}}$,
the matrix ${\tt btph}$ is block diagonalizable by $\left(C_{I_{3}}\otimes
I_{I_{1}}\right)$ and we have

$\left(C_{I_{3}}\otimes I_{I_{1}}\right){\tt
btph}\left(\mathcal{A}\right)\left(C_{I_{3}}^{*}\otimes I_{I_{2}}\right)={\tt
bdiag}\left(\hat{\mathcal{A}}\right),$ (16)

where $\hat{\mathcal{A}}=\mathcal{A}\times_{3}M$ with $M$ is defined by (7)
and

${\tt
bdiag}\left(\hat{\mathcal{A}}\right)=\begin{pmatrix}\hat{\mathcal{A}}^{(1)}&&&\\\
&\hat{\mathcal{A}}^{(2)}&&\\\ &&\ddots&\\\
&&&\hat{\mathcal{A}}^{(I_{3})}\end{pmatrix}\in\mathbb{K}^{I_{1}I_{3}\times
I_{2}I_{3}}.$ (17)

From this last result, we can define the c-product between two third-order
tensors $\mathcal{A}$ and $\mathcal{B}$ of appropriate sizes with $I_{3}$
frontal slices. Before giving the generalized version of the high-order
c-product, we give some definitions and notations. First, we will call a
tensor in the scalar space ($\mathbb{K}_{I_{3}\times\dots\times I_{N}}$) a
scalar-tensor, which will replace the notion of tubes in the case of third-
order tensors. For an $N^{th}$-order tensor
$\mathcal{X}\in\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}\times\dots\times I_{N}}$,
we will define the operator $\tt Vec$ that transforms the tensor into a matrix
as follows

${\tt
Vec(\mathcal{X})}=[X^{1},\,X^{2},\,...,\,X^{n_{3}},\,X^{n_{3}+1},\,...,\,X^{P}]$
(18)

where $P=I_{3}I_{4}\dots I_{N}$ and
$X^{p}=\mathcal{X}(:,:,k_{3},k_{4},\dots,k_{N})\in\mathbb{K}^{I_{1}\times
I_{2}}$ with
$p=k_{3}+\displaystyle\sum_{i=4}^{N}\left(\dfrac{\left(k_{i}-1\right)P}{\prod_{s=i}^{N}I_{i}}\right).$
This generalizes the notion of frontal slices in the case of third-order
tensors. Notice that in the case of scalar-tensors
$\mathcal{X}\in\mathbb{K}^{1\times 1}_{I_{3}\times\dots\times I_{N}}$, the
$X^{p}$’s are scalars and ${\tt Vec(\mathcal{X})}$ is a row vector.
In [2], the authors defined the block-Toeplitz-plus-Hankel matrix for third-
order tensor and here we will define the block-Toelitz-plus-Hankel matrix for
an $N^{th}$-order tensor for $N\geq 3$ by using the block Toeplitz plus Hankel
for an $(N-1)^{th}$-order tensor. To explain this we will present this
procedure only for a fourth-order tensor
$\mathcal{X}\in\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}\times I_{4}}$. We define
the matrix ${\tt btph}\in\mathbb{R}^{n_{1}P\times n_{2}P}$ as the matrix
block-Toeplitz-plus-Hankel where each block of the matrices Toeplitz and
Hankel is a block-Toeplitz-plus-Hankel matrix of third-order tensors,
respectively. Then the block-Toeplitz-plus-Hankel matrix of $\mathcal{X}$ is
defined as follows

$\displaystyle{\tt btph}(\mathcal{X})$ $\displaystyle=$
$\displaystyle\begin{pmatrix}{\tt btph}\left(\mathcal{X}(:,:,:,1)\right)&{\tt
btph}\left(\mathcal{X}(:,:,:,2)\right)&\dots&\dots&{\tt
btph}\left(\mathcal{X}(:,:,:,I_{4})\right)\\\ {\tt
btph}\left(\mathcal{X}(:,:,:,2)\right)&{\tt
btph}\left(\mathcal{X}(:,:,:,1)\right)&\dots&\dots&{\tt
btph}\left(\mathcal{X}(:,:,:,I_{4}-1)\right)\\\
\vdots&\vdots&\vdots&\vdots&\vdots\\\ {\tt
btph}\left(\mathcal{X}(:,:,:,I_{4})\right)&{\tt
btph}\left(\mathcal{X}(:,:,:,I_{4}-1)\right)&\dots&\dots&{\tt
btph}\left(\mathcal{X}(:,:,:,1)\right)\end{pmatrix}$ $\displaystyle+$
$\displaystyle\begin{pmatrix}{\tt
btph}\left(\mathcal{X}(:,:,:,2)\right)&\dots&{\tt
btph}\left(\mathcal{X}(:,:,:,I_{4})\right)&\textbf{0}\\\
\vdots&\reflectbox{$\vdots$}&\reflectbox{$\vdots$}&{\tt
btph}\left(\mathcal{X}(:,:,:,I_{4})\right)\\\ {\tt
btph}\left(\mathcal{X}(:,:,:,I_{4})\right)&\textbf{0}&\reflectbox{$\vdots$}&\vdots\\\
\textbf{0}&{\tt btph}\left(\mathcal{X}(:,:,:,I_{4})\right)&\dots&{\tt
btph}\left(\mathcal{X}(:,:,:,2)\right)\end{pmatrix}.$

${\tt ten}({\tt btph}(\mathcal{X}))=\mathcal{X}$, this operator allows to
reconstruct the original tensor from its associate ${\tt btph}$ matrix. The
block diagonal matrix of the tensor $\mathcal{X}$ is given by

${\tt
bdiag}(\mathcal{X})=\begin{pmatrix}X^{1}&\textbf{0}&\dots&\dots&\textbf{0}\\\
\textbf{0}&X^{2}&\textbf{0}&\dots&\textbf{0}\\\
\vdots&\ddots&\ddots&\ddots&\vdots\\\ \vdots&\ddots&\ddots&\ddots&\vdots\\\
\textbf{0}&\dots&\dots&\dots&X^{P}\end{pmatrix}\in\mathbb{R}^{n_{1}P\times
n_{2}P},$ (19)

where the $X^{i}$ ’s matrices of size $I_{1}\times I_{2}$ are given in (18) as
the representative matrices of $\mathcal{X}$. Notice that in the case of
scalar-tensors, the blocks $X^{i}$ are scalars in $\mathbb{C}$ and in this
case, the block diagonal matrix (19) is just a diagonal matrix and ${\tt
bdiag}(\mathcal{X})={\tt diag}(\mathcal{X})$.

###  3.1 Generalized cosine product c-product

In this subsection, we will introduce a generalized version of the c-product
for high-order tensors. To this end, we first need some theoretical results.
First, remind that in the third-order case, the Toeplitz-plus-Hankel matrix of
a tube $\textbf{a}\in\mathbb{K}_{I_{3}}$ is diagonalizable using the DCT
matrix of order $I_{3}$, i.e.,

$C_{I_{3}}{\tt tph}(\textbf{a})C_{I_{3}}^{*}={\tt
diag}(d);\;where\;d=W_{I_{3}}^{-1}C_{I_{3}}\left(I+Z_{I_{3}}\right){\tt
vec}(\textbf{a}),$ (20)

with $W_{I_{3}}$ and $Z_{I_{3}}$ as defined earlier and ${\tt
vec}(\textbf{a})$ is the vector of size $I_{3}$ whose elements are the
coefficients of the tube a . This result is extended to the high-order case by
considering scalar-tensors instead of tubes.

######  Theorem 6.

Let $\mathcal{X}\in\mathbb{K}_{I_{3}\times\dots\times I_{N}}$, a scalar-
tensor, then its block-Toeplitz-plus-Hankel matrix is diagonalizabl, and we
have

$\left(C_{I_{N}}\otimes C_{I_{N-1}}\otimes\dots\otimes C_{I_{3}}\right){\tt
btph}\left(\mathcal{X}\right)\left(C_{I_{N}}^{*}\otimes
C_{I_{N-1}}^{*}\otimes\dots\otimes C_{I_{3}}^{*}\right)={\tt
diag}\left(\hat{\mathcal{X}}\right)$ (21)

where
$\hat{\mathcal{X}}=\mathcal{X}\times_{3}M_{I_{3}}\times_{4}M_{I_{4}}\times_{5}\dots\times_{N}M_{I_{N}}$
and $C_{I_{N}}$ are DCT matrices.

###### Proof.

For simplicity we consider only the case of fourth-order scalar-tensors
$\mathcal{X}$ in $\mathbb{K}^{1\times 1}_{I_{3}\times I_{4}}$. Using the fact
that the matrix block-Toeplitz-plus-Hankel for a third-order tensor is block-
diagonalizable using the discrete cosine matrix, we get

$\displaystyle\left(I_{I_{4}}\otimes C_{I_{3}}\right){\tt
btph}\left(\mathcal{X}\right)\left(I_{I_{4}}\otimes C_{I_{3}}^{*}\right)={\tt
btph}\left(\mathcal{X}^{\prime}\right),$ (22)

with $\mathcal{X}^{\prime}$ is the tensor such that
$\mathcal{X}^{\prime}(:,:,:,k_{4})=C_{I_{3}}{\tt
btph}(\mathcal{X}(:,:,:,k_{4}))C_{I_{3}}^{*}$ for $k_{4}=1,\dots,\,I_{4}$.
Hence $\mathcal{X}^{\prime}=\mathcal{X}\times_{3}M_{I_{3}}$ and

$\displaystyle{\tt btph}\left(\mathcal{X}^{\prime}\right)$ $\displaystyle=$
$\displaystyle\begin{pmatrix}diag(\mathcal{X}^{\prime}_{1})&\dots&\dots&diag(\mathcal{X}^{\prime
I_{4}})\\\
diag(\mathcal{X}^{\prime}_{2})&diag(\mathcal{X}^{\prime}_{1})&\dots&diag(\mathcal{X}^{\prime}_{I_{4}-1})\\\
\vdots&\ddots&\ddots&\vdots\\\
diag(\mathcal{X}^{\prime}_{I_{4}})&\dots&\dots&diag(\mathcal{X}^{\prime}_{1})\end{pmatrix}$
$\displaystyle+$
$\displaystyle\begin{pmatrix}diag(\mathcal{X}^{\prime}_{2})&\dots&diag(\mathcal{X}^{\prime}_{I_{4}})&\textbf{0}\\\
\vdots&\reflectbox{$\ddots$}&\reflectbox{$\ddots$}&diag(\mathcal{X}^{\prime}_{I_{4}})\\\
diag(\mathcal{X}^{\prime}_{I_{4}})&\textbf{0}&\reflectbox{$\ddots$}&\vdots\\\
\textbf{0}&diag(\mathcal{X}^{\prime}_{I_{4}})&\dots&diag(\mathcal{X}^{\prime}_{2})\end{pmatrix}$

where $\mathcal{X}^{\prime}_{i}=\mathcal{X}^{\prime}(:,:,:,i)$ for
$1=1,2,\dots,I_{4}$. Using a similar result as (16), we get

$\left(C_{I_{4}}\otimes I_{I_{3}}\right){\tt
btph}\left(\mathcal{X}^{\prime}\right)\left(C_{I_{4}}^{*}\otimes
I_{I_{3}}\right)={\tt diag}\left(d\right)\in\mathbb{K}^{I_{3}I_{4}\times
I_{3}I_{4}}$ (23)

Then, from (22) and (23) we obtain

$\displaystyle\left(C_{I_{4}}\otimes I_{I_{3}}\right)\left(I_{I_{4}}\otimes
C_{I_{3}}\right){\tt btph}\left(\mathcal{X}\right)\left(I_{I_{4}}\otimes
C_{I_{3}}^{*}\right)\left(C_{I_{4}}^{*}\otimes I_{I_{3}}\right)={\tt
diag}\left(d\right).$ (24)

Therefore, we can deduce that ${\tt diag}(d)={\tt
bdiag}(\hat{\mathcal{X}})={\tt diag}\left(\hat{\mathcal{X}}\right)$, with
$\hat{\mathcal{X}}=\mathcal{X}\times_{3}M_{I_{3}}\times M_{I_{4}}$ is a
scalar-tensor, with
$M_{I_{i}}=W_{I_{i}}^{-1}C_{I_{i}}\left(I_{i}+Z_{I_{i}}\right)$ for $i=3,4$.
Then, from (24), we obtain

$\displaystyle\left(C_{I_{4}}\otimes C_{I_{3}}\right){\tt
btph}\left(\mathcal{X}\right)\left(C_{I_{4}}^{*}\otimes
C_{I_{3}}^{*}\right)={\tt bdiag}\left(\hat{\mathcal{X}}\right).$

The reslut for fifth-order tensor could be obtained from the one for the
fourth-order case and so on, the result for an $N^{th}$-order scalar-tensor,
can be found recursively from the one for an $(N-1)^{th}$-order scalar-tensor.
This is due to the fact that the $\tt btph$ matrix of an $N^{th}$-order tensor
is comoputed by using the $\tt btph$ matrix of an $(N-1)^{th}$-order tensor.
Therefore, for a general scalar-tensor order
$\mathcal{X}\in\mathbb{K}^{1\times 1}_{I_{3}\times\dots\times_{4}I_{N}}$, we
obtain

$\left(C_{I_{N}}\otimes C_{I_{N-1}}\otimes\dots\otimes C_{I_{3}}\right){\tt
btph}(\mathcal{X})\left(C_{I_{N}}^{*}\otimes
C_{I_{N-1}}^{*}\otimes\dots\otimes C_{I_{3}}^{*}\right)={\tt
diag}(\hat{\mathcal{X}}),$

where
$\hat{\mathcal{X}}=\mathcal{X}\times_{3}M_{I_{3}}\times_{4}M_{I_{4}}\times_{5}\dots\times
M_{I_{N}}$. ∎

Next, we define the new operator $\mathcal{L}$ as follows.

######  Definition 7.

Let $\mathcal{L}:\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}\times\dots\times
I_{N}}\longrightarrow\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}\times\dots\times
I_{N}}$ be the operator defined by

$\mathcal{L}\left(\mathcal{A}\right)=\hat{\mathcal{A}}=\mathcal{A}\times_{3}M_{I_{3}}\times_{4}M_{I_{4}}\times_{5}\dots\times_{N}M_{I_{N}},$

and its inverse

$\mathcal{A}=\mathcal{L}^{-1}\left(\hat{\mathcal{A}}\right)=\hat{\mathcal{A}}\times_{N}M_{I_{N}}^{-1}\times_{N-1}M_{I_{N-1}}^{-1}\times_{N-2}\dots\times_{3}M_{I_{3}}^{-1},$

where $M_{I_{i}}=W_{I_{i}}^{-1}C_{I_{i}}\left(I_{I_{i}}+Z_{I_{i}}\right)$ with
$W_{I_{i}}={\tt diag}(C_{I_{i}}(:,1))$ and $C_{I_{i}}$ is the matrix of
discrete cosine and $Z_{I_{i}}$, $i=3,4,\dots,N$ was already defined.

In the last theorem we proved that a block-Toeplitz-plus-Hankel matrix of a
scalar-tensor is diagonalizable and in the next theorem we will prove that a
block-Toeplitz-plus-Hankel matrix of an $N^{th}$-order tensor is block
diagonalizable.

######  Theorem 8.

Let $\mathcal{X}\in\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}\times\dots\times
I_{N}}$, then its block Toeplitz-plus-Hankel matrix is block diagonalizable
and

${\tt btph}(\mathcal{X})=\left(C_{I_{N}}^{*}\otimes\dots\otimes
C_{I_{3}}^{*}\otimes I_{I_{1}}\right){\tt
bdiag}\left(\mathcal{L}(\mathcal{X})\right)\left(C_{I_{N}}\otimes\dots\otimes
C_{I_{3}}\otimes I_{I_{2}}\right)$ (25)

###### Proof.

For a tensor $\mathcal{X}\in\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}\times\dots\times I_{K}}$ we have

$\displaystyle{\tt
btph}(\mathcal{X})=\sum_{i_{1}=1}^{I_{1}}\sum_{i_{2}=1}^{I_{2}}{\tt
btph}(\mathcal{X}(i_{1},\,i_{2},\,:,\dots,\,:))\otimes
e_{i_{1}}e_{i_{2}}^{T},$

where $e_{i}$ is the $i$-th canonical vector of size $I_{i}$ for $i=1,2$.
Therefore, using Thorem 6, we get

$\displaystyle{\tt btph}(\mathcal{X})$ $\displaystyle=$
$\displaystyle\sum_{i_{1},i_{2}=1}^{I_{1};I_{2}}\left(C_{I_{N}}^{*}\otimes\dots\otimes
C_{I_{3}}^{*}\right){\tt
bdiag}\left(\mathcal{L}(\mathcal{X}(i_{1},\,i_{2},\,:,\,,\dots,\,:))\right)\left(C_{I_{N}}\otimes\dots\otimes
C_{I_{3}}\right)\otimes e_{i_{1}}e_{i_{2}}^{T}$ $\displaystyle=$
$\displaystyle\left(C_{I_{N}}^{*}\otimes\dots\otimes C_{I_{3}}^{*}\otimes
I_{I_{1}}\right){\tt
bdiag}(\mathcal{L}\left(\mathcal{X}\right))\left(C_{I_{N}}\otimes\dots\otimes
C_{I_{3}}\otimes I_{I_{2}}\right).$

∎

Next, we define the generalized c-product.

######  Definition 9.

Let $\mathcal{A}\in\mathbb{K}^{I_{1}\times l}_{I_{3}\times\dots\times I_{N}}$
and $\mathcal{B}\in\mathbb{K}^{l\times I_{2}}_{I_{3}\times\dots\times I_{N}}$.
The generalized c-product $\mathcal{A}*_{c}\mathcal{B}$ is the tensor of
$\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}\times\dots\times I_{N}}$ defined as
follows

$\mathcal{A}*_{c}\mathcal{B}=ten\left({\tt btph}(\mathcal{A}){\tt
btph}(\mathcal{B})\right).$ (26)

In the next definition, we generalize the face-wice product that was already
defined for the third-order tensors.

######  Definition 10.

Generalized face-wise product
Let $\mathcal{A}\in\mathbb{K}^{I_{1}\times l}_{I_{3}\times\dots\times I_{N}}$
and $\mathcal{B}\in\mathbb{K}^{l\times I_{2}}_{I_{3}\times\dots\times I_{N}}$,
we define the face-wise product between $\mathcal{A}$ and $\mathcal{B}$ by the
tensor $\mathcal{C}=\mathcal{A}\triangle\mathcal{B}\in\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}\times\dots\times I_{N}}$, where the $p^{th}$ representative
matrix of $\mathcal{C}$, is computed by the product of the $p^{th}$
representative matrices of $\mathcal{A}$ and $\mathcal{B}$, respectively.
i.e.,

$C^{p}=A^{p}B^{p},\;\;p=1,\,2,\,\dots,\,P,$ (27)

where the matrices $A^{p}$, $B^{p}$ and $C^{p}$ for $p=1,2,\dots,P$ are the
representative matrices given by (18) of $\mathcal{A}$, $\mathcal{B}$ and
$\mathcal{C}$, respectively.

######  Lemma 1.

The generalized c-product of two $N^{th}$-order tensors
$\mathcal{A}\in\mathbb{K}^{I_{1}\times l}_{I_{3}\times\dots\times I_{N}}$ and
$\mathcal{B}\in\mathbb{K}^{l\times I_{2}}_{I_{3}\times\dots\times I_{N}}$ can
be also computed in the cosine domain by

$\mathcal{A}*_{c}\mathcal{B}=\mathcal{L}^{-1}\left(\mathcal{L}(\mathcal{A})\triangle\mathcal{L}(\mathcal{B})\right).$
(28)

###### Proof.

The generalized c-product of the two tensors
$\mathcal{A}\in\mathbb{K}^{I_{1}\times l}_{I_{3}\times\dots I_{N}}$ and
$\mathcal{B}\in\mathbb{K}^{l\times I_{2}}_{I_{3}\times\dots I_{N}}$, is given
by

$\mathcal{A}*_{c}\mathcal{B}=ten\left({\tt btph}(\mathcal{A}){\tt
btph}(\mathcal{B})\right)$

Therefore, using the notation
${\textbf{C}}_{N}^{\,p}=\left(C_{I_{N}}\otimes\dots\otimes C_{I_{3}}\otimes
I_{p}\right)$ for $p\geq 1$ and $N\geq 3$, we get

$\displaystyle{\tt btph}\left(\mathcal{A}*_{c}\mathcal{B}\right)$
$\displaystyle=$ $\displaystyle{\tt btph}\left(\mathcal{A}\right){\tt
btph}\left(\mathcal{B}\right)$ $\displaystyle=$
$\displaystyle\left(\textbf{C}_{N}^{\,I_{1}}\right)^{*}\textbf{C}_{N}^{\,I_{1}}{\tt
btph}\left(\mathcal{A}\right)\left(\textbf{C}_{N}^{\,l}\right)^{*}\textbf{C}_{N}^{\,l}{\tt
btph}\left(\mathcal{B}\right)\left(\textbf{C}_{N}^{I_{2}}\right)^{*}\textbf{C}_{N}^{\,I_{2}}$
$\displaystyle=$ $\displaystyle\left(\textbf{C}_{N}^{\,I_{1}}\right)^{*}{\tt
bdiag}(\mathcal{L}\left(\mathcal{A}\right)){\tt
bdiag}(\mathcal{L}\left(\mathcal{B}\right))\textbf{C}_{N}^{\,I_{2}}$
$\displaystyle=$ $\displaystyle{\tt
btph}\left(\mathcal{L}^{-1}\left(\mathcal{L}\left(\mathcal{A}\right)\triangle\mathcal{L}\left(\mathcal{A}\right)\right)\right),$

It follows that

$\mathcal{A}*_{c}\mathcal{B}=\mathcal{L}^{-1}\left(\mathcal{L}(\mathcal{A})\triangle\mathcal{L}(\mathcal{B})\right).$

∎

###  3.2 Generalized tensor-tensor product

In this subsection, we define a general tensor-tensor $\mathcal{L}$ product
for high-order tensor. We first define $\mathcal{L}$ as the following operator

$\displaystyle\mathcal{L}:\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}\times\dots\times I_{N}}$ $\displaystyle\longrightarrow$
$\displaystyle\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}\times\dots\times I_{N}}$
$\displaystyle\mathcal{A}$ $\displaystyle\longmapsto$
$\displaystyle\mathcal{L}\left(\mathcal{A}\right)=\mathcal{A}\times_{3}M_{3}\dots\times_{N}M_{N}$

with $M_{i}\in\mathbb{K}^{I_{i}\times I_{i}}$ such that
$M_{i}=\alpha_{i}R_{i}$ for $i=3,\,4,\,\dots,\,N$, where $\alpha_{i}>0$ and
$R_{i}$ is an unitary matrix. The inverse operator of $\mathcal{L}$ is defined
as
$\mathcal{L}^{-1}\left(\mathcal{A}\right)=\mathcal{A}\times_{3}M_{3}^{-1}\times_{4}\dots\times_{N}M_{N}^{-1}$.
We will denote $\alpha=\alpha_{3}\alpha_{4}\dots\alpha_{N}.$ Next, we need the
relation between the norm of a tensor and its norm in the transformed domain
(for example Fourier or cosine) given by

$\left\|\mathcal{A}\right\|_{F}=\dfrac{1}{\sqrt{\alpha}}\left\|\mathcal{L}\left(\mathcal{A}\right)\right\|_{F},$
(29)

and we also have

$\left<\mathcal{A},\,\mathcal{B}\right>=\dfrac{1}{\alpha}\left<\mathcal{L}\left(\mathcal{A}\right),\,\mathcal{L}\left(\mathcal{B}\right)\right>,$
(30)

In the next and for any $N^{th}$-order tensor $\mathcal{X}$, we will denote
the $p^{th}$ representative matrix (18) of
$\mathcal{L}\left(\mathcal{X}\right)$ by $\mathcal{L}(\mathcal{X})^{p}$. Now
we can define the generalized $*_{\mathcal{L}}$-product

######  Definition 11.

Let $\mathcal{L}$ be the operator defined above, then the generalized
$*_{\mathcal{L}}$-product of two $N^{th}$-order tensors
$\mathcal{A}\in\mathbb{K}^{I_{1}\times l}_{I_{3}\times\dots I_{N}}$ and
$\mathcal{B}\in\mathbb{K}^{l\times I_{2}}_{I_{3}\times\dots\times I_{N}}$ is
given by

$\mathcal{A}*_{\mathcal{L}}\mathcal{B}=\mathcal{L}^{-1}\left(\mathcal{L}\left(\mathcal{A}\right)\triangle\mathcal{L}\left(\mathcal{B}\right)\right)\in\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}\times\dots\times I_{N}}$ (31)

The whole steps are summarized in the following algorithm.

Algorithm 2 The $*_{\mathcal{L}}$-product.

1: Inputs: $\mathcal{A}\in\mathbb{K}^{n_{1}\times n}_{n_{3}\times\dots\times
n_{N}}$ and $\mathcal{B}\in\mathbb{K}^{n\times n_{2}}_{n_{3}\times\dots\times
n_{N}}$.

2: Output: $\mathcal{C}=\mathcal{A}*\mathcal{B}\in\mathbb{K}^{n_{1}\times
n_{2}}_{n_{3}\dots\times n_{K}}$ .

3: Compute $\hat{\mathcal{A}}=\mathcal{L}\left(\mathcal{A}\right)$ and
$\hat{\mathcal{B}}=\mathcal{L}\left(\mathcal{B}\right)$.

4: for $i=1,\ldots,P$ do

5: $\hat{\mathcal{C}}^{i}=\hat{\mathcal{A}}^{i}\hat{\mathcal{B}}^{i}$

6: end for

7: $\mathcal{C}=\mathcal{L}^{-1}\left(\hat{\mathcal{C}}\right).$

######  Proposition 1.

Let $\mathcal{X}\in\mathbb{K}^{I_{1}\times l}_{I_{3}\times\dots\times I_{K}}$
and $\mathcal{Y}\in\mathbb{K}^{l\times I_{2}}_{I_{3}\times\dots\times I_{K}}$
two $K^{th}$-order tensors. Then we have

$\left(\mathcal{X}*_{\mathcal{L}}\mathcal{Y}\right)_{i,j}=\sum_{k=1}^{l}\mathcal{X}_{i,k}*_{\mathcal{L}}\mathcal{Y}_{k,j},\;\text{for}\;1\leq
i\leq I_{1}\;\text{and}\;1\leq j\leq I_{2},$ (32)

where
$\mathcal{X}_{i,j}=\mathcal{X}(i,j,:,:,\dots,:)\in\mathbb{K}_{I_{3}\times
I_{4}\times\dots\times I_{K}}$ is the $(i,j)^{th}$ scalar-tensor of
$\mathcal{X}$.

###### Proof.

For or $1\leq i\leq I_{1}$ and $1\leq j\leq I_{2}$, we have

$\displaystyle\left(\mathcal{X}*_{\mathcal{L}}\mathcal{Y}\right)_{i,j}$
$\displaystyle=$
$\displaystyle\mathcal{L}^{-1}\left(\mathcal{L}\left(\mathcal{X}\right)\triangle\mathcal{L}\left(\mathcal{Y}\right)\right)_{i,j}$
$\displaystyle=$
$\displaystyle\mathcal{L}^{-1}\left(\sum_{k=1}^{l}\mathcal{L}\left(\mathcal{X}\right)_{i,k}\triangle\mathcal{L}\left(\mathcal{Y}\right)_{k,j}\right)$
$\displaystyle=$
$\displaystyle\sum_{k=1}^{l}\mathcal{L}^{-1}\left(\mathcal{L}\left(\mathcal{X}_{i,k}\right)\triangle\mathcal{L}\left(\mathcal{Y}_{k,j}\right)\right)$
$\displaystyle=$
$\displaystyle\sum_{k=1}^{l}\mathcal{X}_{i,k}*_{\mathcal{L}}\mathcal{Y}_{k,j},$

which shows the result. ∎

######  Proposition 2.

Let $\mathcal{A}\in\mathbb{K}^{I_{1}\times l}_{I_{3}\times\dots\times I_{K}}$
and $\mathcal{B}\in\mathbb{K}^{l\times I_{2}}_{I_{3}\times\dots\times I_{K}}$
be two $K^{th}$-order tensors, then we can express the
$*_{\mathcal{L}}$-product of $\mathcal{A}$ and $\mathcal{B}$ as

$\mathcal{A}*_{\mathcal{L}}\mathcal{B}=\left[\mathcal{A}*_{\mathcal{L}}\overrightarrow{\mathcal{B}}_{1},\,\mathcal{A}*_{\mathcal{L}}\overrightarrow{\mathcal{B}}_{2},\,\dots,\,\mathcal{A}*_{\mathcal{L}}\overrightarrow{\mathcal{B}}_{I_{2}}\right],$
(33)

where
$\overrightarrow{\mathcal{B}}_{i_{2}}=\mathcal{B}(:,i_{2},:,\dots,:)\in\mathbb{K}^{l}_{I_{3}\times
I_{4}\times\dots\times I_{K}};\;i_{2}=1,\ldots,I_{2}$.

###### Proof.

From the definition of the $*_{\mathcal{L}}$-product of $\mathcal{A}$ and
$\mathcal{B}$ and for $1\leq i_{2}\leq I_{2}$, we get

$\displaystyle\overrightarrow{\left(\mathcal{A}*_{\mathcal{L}}\mathcal{B}\right)}_{i_{2}}$
$\displaystyle=$
$\displaystyle\overrightarrow{\mathcal{L}^{-1}\left(\mathcal{L}\left(\mathcal{A}\right)\triangle\mathcal{L}\left(\mathcal{B}\right)\right)}_{i_{2}}$
(34) $\displaystyle=$
$\displaystyle\mathcal{L}^{-1}\left(\mathcal{L}\left(\mathcal{A}\right)\triangle\overrightarrow{\mathcal{L}\left(\mathcal{B}\right)}_{i_{2}}\right)$
(35) $\displaystyle=$
$\displaystyle\mathcal{L}^{-1}\left(\mathcal{L}\left(\mathcal{A}\right)\triangle\mathcal{L}\left(\overrightarrow{\mathcal{B}}_{i_{2}}\right)\right)$
(36) $\displaystyle=$
$\displaystyle\mathcal{A}*_{\mathcal{L}}\overrightarrow{\mathcal{B}}_{i_{2}},$
(37)

which gives the desired result. ∎

Related to the generalized $*_{\mathcal{L}}$-product , we give the definitions
of the identity, transpose and orthogonal tensors.

######  Definition 12.

(The identity tensor)
The tensor identity tensor $\mathcal{I}\in\mathbb{K}^{I_{1}\times
I_{1}}_{I_{3}\times\dots I_{N}}$ is such that
$\left(\mathcal{L}(\mathcal{I})\right)^{p}=I$ for $p=1,2,\dots,P$ where
$\left(\mathcal{L}(\mathcal{I})\right)^{p}$ is the $p^{th}$ representative
matrix of $\mathcal{L}\left(\mathcal{I}\right)$.

From the previous definition, we can conclude that if
$\mathcal{A}\in\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}\times\dots\times I_{N}}$,
then
$\mathcal{A}*_{\mathcal{L}}\mathcal{I}=\mathcal{I}*_{\mathcal{L}}\mathcal{A}=\mathcal{A}$;
because
$\mathcal{A}*_{\mathcal{L}}\mathcal{I}=\mathcal{L}^{-1}\left(\mathcal{L}(\mathcal{A})\triangle\mathcal{L}(\mathcal{I})\right)=\mathcal{L}^{-1}\left(\mathcal{L}(\mathcal{A})\right),$
and the same for $\mathcal{I}*_{\mathcal{L}}\mathcal{A}.$

######  Definition 13.

(The transpose)
Let $\mathcal{A}\in\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}\times\dots\times
I_{N}}$, the transpose $\mathcal{A}^{T}$ of the tensor $\mathcal{A}$ is such
that
$\mathcal{L}\left(\mathcal{A}^{T}\right)^{p}=\left(\mathcal{L}\left(\mathcal{A}\right)^{p}\right)^{T}$
for $p=1,2,\dots,P$.

This definition ensure the multiplication reversal property for the transpose
under the $*_{\mathcal{L}}$-product, i.e., for $\mathcal{A}$ and $\mathcal{B}$
two $N^{th}$-order tensors of appropriate sizes, we get
$\left(\mathcal{A}*_{\mathcal{L}}\mathcal{B}\right)^{T}=\mathcal{B}^{T}*_{\mathcal{L}}\mathcal{A}^{T}$;
for explanation we have

$\displaystyle\mathcal{L}\left(\left(\mathcal{A}*_{\mathcal{L}}\mathcal{B}\right)^{T}\right)^{p}=\left(\mathcal{L}\left(\mathcal{A}*_{\mathcal{L}}\mathcal{B}\right)^{p}\right)^{T}$
$\displaystyle=$
$\displaystyle\left(\mathcal{L}\left(\mathcal{A}\right)^{p}\mathcal{L}\left(\mathcal{B}\right)^{p}\right)^{T}$
$\displaystyle=$
$\displaystyle\left(\mathcal{L}\left(\mathcal{B}\right)^{p}\right)^{T}\left(\mathcal{L}\left(\mathcal{A}\right)^{p}\right)^{T}$
$\displaystyle=$
$\displaystyle\left(\mathcal{L}\left(\mathcal{B}^{T}\right)\right)^{p}\left(\mathcal{L}\left(\mathcal{A}^{T}\right)\right)^{p}$
$\displaystyle=$
$\displaystyle\mathcal{L}\left(\mathcal{B}^{T}*_{\mathcal{L}}\mathcal{A}^{T}\right)^{p}$

######  Definition 14.

(Orthogonal tensor)
The tensor $\mathcal{Q}\in\mathbb{K}^{I_{1}\times
I_{1}}_{I_{3}\times\dots\times I_{N}}$ is orthogonal under the
$*_{\mathcal{L}}$-product, iff
$\mathcal{Q}*_{\mathcal{L}}\mathcal{Q}^{T}=\mathcal{Q}^{T}*_{\mathcal{L}}\mathcal{Q}=\mathcal{I}$
which means that for each $i\in\\{1,\,2,\,\dots,\,P\\}$,
$\mathcal{L}\left(\mathcal{Q}\right)^{i}$ is an orthogonal matrix.

Notice that if $\mathcal{Q}\in\mathbb{K}^{I_{1}\times
I_{1}}_{I_{3}\times\dots\times I_{N}}$ is orthogonal, then for an
$N^{th}$-order tensor $\mathcal{A}$ of an appropriate size, we have

$\displaystyle\left\|\mathcal{A}*_{\mathcal{L}}\mathcal{Q}\right\|_{F}^{2}=\dfrac{1}{c}\left\|\mathcal{L}\left(\mathcal{A}\right)\triangle\mathcal{L}\left(\mathcal{Q}\right)\right\|_{F}^{2}$
$\displaystyle=$
$\displaystyle\dfrac{1}{c}\sum_{i=1}^{P}\left\|\mathcal{L}\left(\mathcal{A}\right)^{i}\mathcal{L}\left(\mathcal{Q}\right)^{i}\right\|_{F}^{2}$
(38) $\displaystyle=$
$\displaystyle\dfrac{1}{c}\sum_{i=1}^{P}\left\|\mathcal{L}\left(\mathcal{A}\right)^{i}\right\|_{F}^{2}$
(39) $\displaystyle=$ $\displaystyle\left\|\mathcal{A}\right\|_{F}^{2}$ (40)

######  Definition 15.

f-diagonal tensor
An $N^{th}$-order tensor $\mathcal{X}$ is f-diagonal, if each
$\mathcal{L}(\mathcal{X})^{p}$ is a diagonal matrix for all $p$ in
$\\{1,2,\dots,P\\}$.

######  Theorem 16.

The set $\left(\mathbb{K}_{I_{3}\times I_{4}\times\dots\times
I_{N}},+,*_{\mathcal{L}}\right)$ is a commutaive ring.

###### Proof.

It is easy to prove that $\left(\mathbb{K}_{I_{3}\times I_{4}\times\dots\times
I_{N}},+\right)$ is an abelian group with
$\textbf{0}\in\mathbb{K}_{I_{3}\times I_{4}\times\dots\times I_{N}}$ as a
neutral element. On the other hand, we have

  * •

Let $\mathcal{A},\,\mathcal{B}$ and $\mathcal{C}$ are scalar-tensors in
$\mathbb{K}_{I_{3}\times I_{4}\times\dots\times I_{N}}$, then we will get

$\displaystyle\mathcal{A}*_{\mathcal{L}}\left(\mathcal{B}*_{\mathcal{L}}\mathcal{C}\right)$
$\displaystyle=$
$\displaystyle\mathcal{L}^{-1}\left(\mathcal{L}\left(\mathcal{A}\right)\triangle\mathcal{L}\left(\mathcal{L}^{-1}\left(\mathcal{L}\left(\mathcal{B}\right)\triangle\mathcal{L}\left(\mathcal{C}\right)\right)\right)\right)$
$\displaystyle=$
$\displaystyle\mathcal{L}^{-1}\left(\mathcal{L}\left(\mathcal{A}\right)\triangle\mathcal{L}\left(\mathcal{B}\right)\triangle\mathcal{L}\left(\mathcal{C}\right)\right)$
$\displaystyle=$
$\displaystyle\mathcal{L}^{-1}\left(\mathcal{L}\left(\mathcal{L}^{-1}\left(\mathcal{L}\left(\mathcal{A}\right)\triangle\mathcal{L}\left(\mathcal{B}\right)\right)\right)\triangle\mathcal{L}\left(\textbf{c}\right)\right)$
$\displaystyle=$
$\displaystyle\left(\mathcal{A}*_{\mathcal{L}}\mathcal{B}\right)*_{\mathcal{L}}\mathcal{C}.$

  * •

We also have

$\displaystyle\mathcal{A}*_{\mathcal{L}}\left(\mathcal{B}+\mathcal{C}\right)$
$\displaystyle=$
$\displaystyle\mathcal{L}^{-1}\left(\mathcal{L}\left(\mathcal{A}\right)\triangle\mathcal{L}\left(\mathcal{B}+\mathcal{C}\right)\right)$
$\displaystyle=$
$\displaystyle\mathcal{L}^{-1}\left(\mathcal{L}\left(\mathcal{A}\right)\triangle\left(\mathcal{L}\left(\mathcal{B}\right)+\mathcal{L}\left(\mathcal{C}\right)\right)\right)$
$\displaystyle=$
$\displaystyle\mathcal{L}^{-1}\left(\mathcal{L}\left(\mathcal{A}\right)\triangle\mathcal{L}\left(\mathcal{B}\right)+\mathcal{L}\left(\mathcal{A}\right)\triangle\mathcal{L}\left(\mathcal{C}\right)\right)$
$\displaystyle=$
$\displaystyle\mathcal{L}^{-1}\left(\mathcal{L}\left(\mathcal{A}\right)\triangle\mathcal{L}\left(\mathcal{B}\right)\right)+\mathcal{L}^{-1}\left(\mathcal{L}\left(\mathcal{A}\right)\triangle\mathcal{L}\left(\mathcal{C}\right)\right)$
$\displaystyle=$
$\displaystyle\mathcal{A}*_{\mathcal{L}}\mathcal{B}+\mathcal{A}*_{\mathcal{L}}\mathcal{C}.$

  * •

Finally, the commutativity is also satisfied by

$\displaystyle\mathcal{A}*_{\mathcal{L}}\mathcal{B}=\mathcal{L}^{-1}\left(\mathcal{L}\left(\mathcal{A}\right)\triangle\mathcal{L}\left(\mathcal{B}\right)\right)$
$\displaystyle=$
$\displaystyle\mathcal{L}^{-1}\left(\mathcal{L}\left(\mathcal{B}\right)\triangle\mathcal{L}\left(\mathcal{A}\right)\right)$
$\displaystyle=$ $\displaystyle\mathcal{B}*_{\mathcal{L}}\mathcal{A}.$

∎

######  Theorem 17.

(The tensor $*_{\mathcal{L}}$-SVD)
Let $\mathcal{A}\in\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}\times\dots\times
I_{N}}$ be an $N^{th}$-order tensor, then $\mathcal{A}$ can be decomposed as

$\mathcal{A}=\mathcal{U}*_{\mathcal{L}}\mathcal{S}*_{\mathcal{L}}\mathcal{V}^{T}=\sum_{k=1}^{r}\mathcal{U}(:,k,:,\dots,:)*_{\mathcal{L}}\mathcal{S}(k,k,:,\dots,:)*_{\mathcal{L}}\mathcal{V}(:,k,:,\dots,:)^{T},$
(41)

where $\mathcal{U}\in\mathbb{K}^{I_{1}\times I_{1}}_{I_{3}\times\dots\times
I_{N}}$ and $\mathcal{V}\in\mathbb{K}^{I_{2}\times
I_{2}}_{I_{3}\times\dots\times I_{N}}$ are orthogonal tensors,
$\mathcal{S}\in\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}\times\dots\times I_{N}}$
is an f-diagonal and $r$ is the tubal rank of $\mathcal{A}$ which will be
defined in the next.

###### Proof.

For a $K^{th}$-order tensor $\mathcal{A}$ and $p=1,\,2,\,\dots,\,P$, consider
SVD decomposition of the matrices

$\mathcal{L}\left(\mathcal{A}\right)^{p}=U_{\mathcal{L}}^{p}S_{\mathcal{L}}^{p}\left(V_{\mathcal{L}}^{p}\right)^{T}.$

Then we obtain

$\mathcal{L}\left(\mathcal{A}\right)=\mathcal{U}_{\mathcal{L}}\triangle\mathcal{S}_{\mathcal{L}}\triangle\left(\mathcal{V}_{\mathcal{L}}\right)^{T},$

with $\mathcal{L}\left(\mathcal{U}\right)=\mathcal{U}_{\mathcal{L}}$,
$\mathcal{L}\left(\mathcal{S}\right)=\mathcal{S}_{\mathcal{L}}$ and
$\mathcal{L}\left(\mathcal{V}\right)=\mathcal{V}_{\mathcal{L}}$, where those
tensors are well defined since the operator $\mathcal{L}$ is inverstible.
Therefore

$\mathcal{L}\left(\mathcal{A}\right)=\mathcal{L}\left(\mathcal{U}\right)\triangle\mathcal{L}\left(\mathcal{S}\right)\triangle\mathcal{L}\left(\mathcal{V}\right)^{T}\Longleftrightarrow\mathcal{A}=\mathcal{L}^{-1}\left(\mathcal{L}\left(\mathcal{U}\right)\triangle\mathcal{L}\left(\mathcal{S}\right)\triangle\mathcal{L}\left(\mathcal{V}\right)^{T}\right)$

which gives

$\mathcal{A}=\mathcal{U}*_{\mathcal{L}}\mathcal{S}*_{\mathcal{L}}\mathcal{V}^{T}.$

∎

The following algorithm summarises the different steps for computing the
tensor $*_{\mathcal{L}}$-SVD of an $N^{th}$-order tensor.

Algorithm 3 The $*_{\mathcal{L}}$-svd.

1: Inputs: $\mathcal{A}\in\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}\times\dots\times I_{N}}$.

2: Output: $\mathcal{U}\in\mathbb{K}^{I_{1}\times
I_{1}}_{I_{3}\times\dots\times I_{N}}$, $\mathcal{S}\in\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}\times\dots\times I_{N}}$, $\mathcal{V}\in\mathbb{K}^{I_{2}\times
I_{2}}_{I_{3}\times\dots\times I_{N}}$.

3: Compute $\hat{\mathcal{A}}=\mathcal{L}\left(\mathcal{A}\right)$.

4: for $i=1,\ldots,P$ do

5:
$\left[\hat{\mathcal{U}}^{i},\hat{\mathcal{S}}^{i},\hat{\mathcal{V}}^{i}\right]=svd(\hat{\mathcal{A}}^{i})$

6: end for

7: $\mathcal{U}=\mathcal{L}^{-1}\left(\hat{\mathcal{U}}\right)$,
$\mathcal{S}=\mathcal{L}^{-1}\left(\hat{\mathcal{S}}\right)$ and
$\mathcal{V}=\mathcal{L}^{-1}\left(\hat{\mathcal{V}}\right)$

######  Corollary 1.

Let $\mathcal{A}$ be an $N^{th}$-order tensor in $\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}\times\dots\times I_{N}}$, and let
$\mathcal{S}_{i}=\mathcal{S}(i,i,:,\dots,:)$, where $\mathcal{S}$ is given by
(41). Then

$\left\|\mathcal{A}\right\|_{F}^{2}=\left\|\mathcal{S}\right\|_{F}^{2}=\sum_{i=1}^{\min(I_{1},\,I_{2})}\left\|\mathcal{S}_{i}\right\|_{F}^{2}.$
(42)

Furthermore

$\left\|\mathcal{S}_{1}\right\|_{F}\geq\left\|\mathcal{S}_{2}\right\|_{F}\geq\dots\geq
0.$ (43)

###### Proof.

Since the tensors $\mathcal{U}$ and $\mathcal{V}$ given from the
$*_{\mathcal{L}}$-svd are orthogonal, then

$\displaystyle\left\|\mathcal{A}\right\|_{F}^{2}=\left\|\mathcal{U}*_{\mathcal{L}}\mathcal{S}*_{\mathcal{L}}\mathcal{V}^{T}\right\|_{F}^{2}$
$\displaystyle=$ $\displaystyle\left\|\mathcal{S}\right\|_{F}^{2}$
$\displaystyle=$
$\displaystyle\dfrac{1}{\alpha}\sum_{p=1}^{P}\left\|\mathcal{L}\left(\mathcal{S}\right)^{p}\right\|_{F}^{2}$
$\displaystyle=$
$\displaystyle\dfrac{1}{\alpha}\sum_{p=1}^{P}\sum_{i=1}^{\min(I_{1},\,I_{2})}\left(\mathcal{L}\left(\mathcal{S}\right)^{p}(i,i)\right)^{2}$
$\displaystyle=$
$\displaystyle\sum_{i=1}^{\min(I_{1},\,I_{2})}\left\|\mathcal{S}_{i}\right\|_{F}^{2}.$

On the other hand,

$\displaystyle\left\|\mathcal{S}_{i}\right\|_{F}^{2}=\dfrac{1}{\alpha}\left\|\mathcal{L}\left(\mathcal{S}\right)_{i}\right\|_{F}^{2}=\dfrac{1}{\alpha}\sum_{p=1}^{P}\left(\mathcal{L}\left(\mathcal{S}\right)^{p}(i,i)\right)^{2}\geq\dfrac{1}{\alpha}\sum_{p=1}^{P}\left(\mathcal{L}\left(\mathcal{S}\right)^{p}(i+1,i+1)\right)^{2}=\left\|\mathcal{S}_{i+1}\right\|_{F}^{2}.$

∎

Next, we give different definitions of a rank of a high-order tensor.

######  Definition 18.

(The $*_{\mathcal{L}}$-tubal rank)
Let $\mathcal{A}\in\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}\times\dots\times
I_{N}}$, then the tensor $*_{\mathcal{L}}$-tubal rank is defined as

${\tt rank}_{t}(\mathcal{A})={\tt
card}\\{i/\;\mathcal{S}(i,i,:,\dots,:)\neq\textbf{0}\\},$ (44)

where $\mathcal{S}$ is the $f$-diagonal tensor given from the
$*_{\mathcal{L}}$-svd of $\mathcal{A}$ (41).

######  Lemma 2.

For a $K^{th}$-order tensor, the $*_{\mathcal{L}}$-tubal rank can be written
as

${\tt rank}_{t}\left(\mathcal{A}\right)={\tt card}\left\\{i/\exists
p\in\\{1,\,2,\,\dots,\,P\\};\;\mathcal{L}\left(\mathcal{S}\right)^{p}(i,i)\neq
0\right\\}.$ (45)

###### Proof.

The key idea for proving the above result is the following equivalence,

$\mathcal{S}\left(i,i,:,\dots,:\right)=\textbf{0}\Longleftrightarrow\forall
p\in\\{1,2,\dots,P\\};\;\mathcal{L}\left(\mathcal{S}\right)^{p}(i,i)=0,\,i=1,2,\dots,\min(I_{1},I_{2}),$

which is equivalent to

$\mathcal{S}\left(i,i,:,\dots,:\right)\neq\textbf{0}\Longleftrightarrow\exists
p\in\\{1,2,\dots,P\\};\;\mathcal{L}\left(\mathcal{S}\right)^{p}(i,i)\neq
0,\;,\,i=1,2,\dots,\min(I_{1},I_{2}).$

Therefore,

$\displaystyle\mathcal{S}\left(i,i,:,\dots,:\right)\neq\textbf{0}\Longleftrightarrow\underset{1\leq
p\leq P}{\max}\,\mathcal{L}\left(\mathcal{S}\right)^{p}(i,i)\neq
0,\,i=1,2,\dots,\min(I_{1},I_{2}).$ (46)

∎

######  Definition 19.

(Multirank and average rank)
For an $N^{th}$-order tensor $\mathcal{A}$ of size $I_{1}\times
I_{2}\times\dots\times I_{N}$, its multirank under the
$*_{\mathcal{L}}$-product is defined as the vector $\rho$ of size $P$, where
its $i^{th}$ element is the rank of $\mathcal{L}\left(\mathcal{A}\right)^{i}$,
i.e.,

$\rho_{i}={\tt rank}(\mathcal{L}\left(\mathcal{A}\right)^{i}).$ (47)

The average rank of $\mathcal{A}$ is defined as the mean of the vector $\rho$,
i.e.,

${\tt
rank}_{a}\left(\mathcal{A}\right)=\dfrac{\displaystyle\sum_{i=1}^{P}\rho_{i}}{P}.$
(48)

######  Remark 1.

We notice that the average rank of an $N^{th}$-order tensor is defined as the
rank of the block-diagonal matrix of $\mathcal{L}\left(\mathcal{A}\right)$
divided by $P$,

${\tt rank}_{a}(\mathcal{A})=\dfrac{{\tt
bdiag}(\mathcal{L}\left(\mathcal{A}\right))}{P}.$ (49)

Next, we give a generalized version of the well known Eckart Young using the
$*_{\mathcal{L}}$-product .

######  Theorem 20.

(Eckart Young)
Let $\mathcal{A}\in\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}\times\dots\times
I_{N}}$ and
$\mathcal{A}_{k}=\mathcal{U}(:,1:k,:,\dots,:)*_{\mathcal{L}}\mathcal{S}(1:k,1:k,:,\dots,:)*_{\mathcal{L}}\mathcal{V}(1:k,:,\dots,:)^{T}$,
where $k$ is a nonzero positive integer. Let $\mathcal{C}$ be the set befined
by

$\mathcal{C}=\left\\{\mathcal{X}*_{\mathcal{L}}\mathcal{Y}\;/\mathcal{X}\in\mathbb{K}^{I_{1}\times
k}_{I_{3}\times\dots\times I_{N}},\;\mathcal{Y}\in\mathbb{K}^{k\times
I_{2}}_{I_{3}\times\dots\times I_{N}}\right\\}.$

Then the tensor $\mathcal{A}_{k}$ solves the minimisation problem

$\displaystyle\min_{\mathcal{X}\in\mathtt{C}}\|\mathcal{A}-\mathcal{X}\|_{F},$

and the error-norm is given by

$\left\|\mathcal{A}-\mathcal{A}_{k}\right\|_{F}^{2}=\sum_{i=k+1}^{r}\left\|\mathcal{S}_{i}\right\|_{F}^{2}.$

###### Proof.

Let $\mathcal{B}=\mathcal{X}*_{\mathcal{L}}\mathcal{Y}\in\mathcal{C}$, then we
have

$\left\|\mathcal{A}-\mathcal{B}\right\|_{F}^{2}=\dfrac{1}{\alpha}\sum_{i=1}^{P}\left\|\mathcal{L}\left(\mathcal{A}\right)^{i}-\mathcal{L}\left(\mathcal{B}\right)^{i}\right\|_{F}^{2}.$

Now, since
$\mathcal{L}\left(\mathcal{B}\right)^{i}=\mathcal{L}\left(\mathcal{X}\right)^{i}\mathcal{L}\left(\mathcal{Y}\right)^{i}$
and by using the matrix Eckart Young theorem we obtain the result showing that
the $k$-best approximation of the matrix
$\mathcal{L}\left(\mathcal{B}\right)^{i}$ is given by
$\mathcal{L}\left(\mathcal{U}\right)^{i}(:,1:k)\mathcal{L}\left(\mathcal{S}\right)^{i}(1:k,1:k)\mathcal{L}\left(\mathcal{V}\right)^{i}(1:k,:)^{T}$,
where $\mathcal{U}$, $\mathcal{S}$ and $\mathcal{V}$ are given by (41). ∎

######  Proposition 3.

The tensor spectral norm of an $N^{th}$-order tensor $\mathcal{A}$ satisfies
the following equation

$\left\|\mathcal{A}\right\|=\left\|{\tt
bdiag}\left(\mathcal{L}\left(\mathcal{A}\right)\right)\right\|$ (50)

###### Proof.

The spectral norm of a tensor $\mathcal{A}$, is defined by

$\displaystyle\left\|\mathcal{A}\right\|=\underset{\underset{\left\|\mathcal{V}\right\|_{F}=1}{\mathcal{V}\in\mathbb{K}^{I_{2}\times
1}_{I_{3}\times\dots\times
I_{N}}}}{sup}\left\|\mathcal{A}*_{\mathcal{L}}\mathcal{V}\right\|_{F}$
$\displaystyle=$
$\displaystyle\dfrac{1}{\sqrt{\alpha}}\underset{\underset{\left\|\mathcal{V}\right\|_{F}=1}{\mathcal{V}\in\mathbb{K}^{I_{2}\times
1}_{I_{3}\times\dots\times
I_{N}}}}{sup}\left\|\mathcal{L}(\mathcal{A})\triangle\mathcal{L}(\mathcal{V})\right\|_{F}$
(51) $\displaystyle=$
$\displaystyle\dfrac{1}{\sqrt{\alpha}}\underset{\underset{\left\|\mathcal{V}\right\|_{F}=1}{\mathcal{V}\in\mathbb{K}^{I_{2}\times
1}_{I_{3}\times\dots\times I_{N}}}}{sup}\left\|{\tt
bdiag}\left(\mathcal{L}(\mathcal{A})\right){\tt
bdiag}\left(\mathcal{L}(\mathcal{V})\right)\right\|_{F}$ (52) $\displaystyle=$
$\displaystyle\left\|{\tt bdiag}\left(\mathcal{L}(\mathcal{A})\right)\right\|$
(53)

∎

######  Definition 21.

Tensor nuclear norm
The tensor nuclear norm of a tensor $\mathcal{A}\in\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}\times\dots\times I_{N}}$ is defined as its dual norm, i.e.,

$\left\|\mathcal{A}\right\|_{*}=\underset{\left\|\mathcal{B}\right\|\leq
1}{sup}\,\left|\left<\mathcal{A},\mathcal{B}\right>\right|$ (54)

######  Theorem 22.

The nuclear norm of a tensor $\mathcal{A}\in\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}\times\dots\times I_{N}}$ verifies the relation

$\left\|\mathcal{A}\right\|_{*}=\dfrac{1}{\alpha}\left\|{\tt
bdiag}(\mathcal{L}\left(\mathcal{A}\right))\right\|_{*}=\dfrac{1}{\alpha}\left\|{\tt
bdiag}(\mathcal{L}\left(\mathcal{S}\right))\right\|_{*}.$ (55)

###### Proof.

Starting by the definition of the nuclear norm in (54), we will get

$\displaystyle\left\|\mathcal{A}\right\|_{*}=\underset{\left\|\mathcal{B}\right\|\leq
1}{\sup}\,\left|\left<\mathcal{A},\mathcal{B}\right>\right|$ $\displaystyle=$
$\displaystyle\dfrac{1}{\alpha}\underset{\left\|\mathcal{B}\right\|\leq
1}{\sup}\,\left|\left<{\tt
bdiag}\left(\mathcal{L}\left(\mathcal{A}\right)\right),{\tt
bdiag}\left(\mathcal{L}\left(\mathcal{B}\right)\right)\right>\right|$
$\displaystyle=$ $\displaystyle\dfrac{1}{\alpha}\left\|{\tt
bdiag}(\mathcal{L}\left(\mathcal{A}\right))\right\|_{*}.$

Which gives the first equality of (55), and the last equality is trivial since
${\tt bdiag}\left(\mathcal{L}\left(\mathcal{A}\right)\right)={\tt
bdiag}\left(\mathcal{L}\left(\mathcal{U}\right)\right){\tt
bdiag}\left(\mathcal{L}\left(\mathcal{S}\right)\right){\tt
bdiag}\left(\mathcal{L}\left(\mathcal{V}^{T}\right)\right)$ where
$\mathcal{U}$, $\mathcal{S}$ and $\mathcal{V}$ are given from the
$*_{\mathcal{L}}$-svd of $\mathcal{A}$. ∎

We can also express the tensor nuclear norm of a $K^{th}$-order tensor as

$\left\|\mathcal{A}\right\|_{*}=\dfrac{1}{\alpha}\sum_{p=1}^{P}\left(\sum_{i=1}^{\min(I_{1},I_{2})}\mathcal{L}\left(\mathcal{S}\right)^{p}(i,i)\right).$
(56)

######  Theorem 23.

The envelope convex of the function tensor average rank on the set
$\mathbb{S}=\left\\{\mathcal{A}\in\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}\times\dots\times I_{N}}/\,\left\|\mathcal{A}\right\|\leq
1\right\\}$ is the tensor nuclear norm.

###### Proof.

Let $\mathcal{A}\in\mathbb{S}$, then $\left\|\mathcal{A}\right\|\leq 1$, which
means that for $i=1,\,2,\,\dots,\,P$, we have

$\left\|\mathcal{L}\left(\mathcal{A}\right)^{i}\right\|\leq 1$

and by using the fact that the tensor average rank is the average of the
tensor multirank, it follows (see [6]) that for each $i=1,\,2,\,\dots,\,P$,
the envelope convex of $\rho_{i}$ is
$\left\|\mathcal{L}\left(\mathcal{A}\right)^{i}\right\|_{*}$.
Consequently, the envelope convex of the function average rank is the tensor
nuclear norm. ∎

Next, we define the tensor singular value thresholding under the
$*_{\mathcal{L}}$-product ($*_{\mathcal{L}}$-svt).

######  Definition 24.

Let $\mathcal{A}\in\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}\times\dots\times
I_{N}}$ an $N^{th}$-order tensor and $\tau>0$, then we call its tensor
singular value thresholding the following tensor

$\mathcal{D}_{\tau}\left(\mathcal{A}\right)=\mathcal{U}*\mathcal{S}_{\tau}*\mathcal{V}^{T}$
(57)

where $\mathcal{U}$, $\mathcal{S}$ and $\mathcal{V}$ are given from the
geberalize t-svd (41), and
$\mathcal{S}_{\tau}=\mathcal{L}^{-1}\left((\mathcal{L}(\mathcal{S})-\tau)_{+}\right)$.

We also have the following important result that links the tensor nuclear norm
and the tensor singular value thresholding.

######  Theorem 25.

Let $\mathcal{A}\in\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}\times\dots\times
I_{N}}$ and $\tau>0$, then we have

$\displaystyle\mathcal{D}_{\tau}\left(\mathcal{A}\right)=\underset{\mathcal{Y}\in\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}\times\dots\times
I_{N}}}{\arg\,\min}\,\tau\left\|\mathcal{Y}\right\|_{*}+\dfrac{1}{2}\left\|\mathcal{Y}-\mathcal{A}\right\|_{F}^{2}$
(58)

###### Proof.

Solving the optimization problem (58) is equivalent to solve the following one

$\underset{\mathcal{Y}\in\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}\times\dots\times
I_{N}}}{\arg\min}\,\dfrac{1}{\alpha}\left\\{\tau\left\|\ {\tt
bdiag}(\mathcal{\mathcal{L}}\left({\mathcal{Y}}\right))\right\|_{*}+\dfrac{1}{2}\left\|{\tt
bdiag}(\mathcal{\mathcal{L}}\left({\mathcal{Y}}\right))-{\tt
bdiag}(\mathcal{\mathcal{L}}\left({\mathcal{A}}\right))\right\|_{F}^{2}\right\\}.$
(59)

The optimization problem given in (59) can be solved by solving $P$
subproblems independently, i.e., for each $i\in\\{1,\,2,\,\dots,\,P\\}$ we
will try to find the solution of

$\underset{\mathcal{\mathcal{L}}\left({\mathcal{Y}}\right)^{i}\in\mathbb{K}^{I_{1}\times
I_{2}}}{\arg\min}\,\left\\{\tau\left\|\mathcal{\mathcal{L}}\left({\mathcal{Y}}\right)^{i}\right\|_{*}+\dfrac{1}{2}\left\|\mathcal{\mathcal{L}}\left({\mathcal{Y}}\right)^{i}-\mathcal{\mathcal{L}}\left({\mathcal{A}}\right)^{i}\right\|_{F}^{2}\right\\}.$
(60)

Using [11], the solution of each subproblem in (60) is
$\mathcal{L}\left(\mathcal{D}_{\tau}\left(\mathcal{A}\right)\right)^{i}$.
Therefore, $\mathcal{D}_{\tau}\left(\mathcal{A}\right)$ solves (58). ∎

In the following algorithm we give the different steps of computing the tensor
$*_{\mathcal{L}}$-svt

Algorithm 4 The $*_{\mathcal{L}}$-svt.

1: Inputs: $\mathcal{A}\in\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}\times\dots\times I_{N}}$ and $\tau>0$.

2: Output:
$\mathcal{D}_{\tau}\left(\mathcal{A}\right)\in\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}\times\dots\times I_{N}}$.

3: Compute: $\hat{\mathcal{A}}=\mathcal{L}\left(\mathcal{A}\right)$

4: for $i=1,\ldots,P$ do

5:
$\left[\hat{\mathcal{U}}^{i},\hat{\mathcal{S}}^{i},\hat{\mathcal{V}}^{i}\right]=svd(\hat{\mathcal{A}}^{i})$

6: $\hat{\mathcal{S}}^{i}=\max(\hat{\mathcal{S}}^{i}-\tau,0)$.

7:
$\hat{\mathcal{A}}^{i}=\hat{\mathcal{U}}^{i}\hat{\mathcal{S}}^{i}(\hat{\mathcal{V}}^{i})^{T}$

8: end for

9: $\mathcal{A}=\mathcal{L}^{-1}(\hat{\mathcal{A}})$.

##  4 Tensor completion using $*_{\mathcal{L}}$-product

Tensor completion is the problem that consists in finding some unknown pixels
of the data from an observed data that contains some known pixels. Many
algorithms have been developed the last years; see [8, 22, 26]. Some of those
methods use regularization techniques such as the total variation
regularization [24, 18]. Those algorithms suffer from the computationally
costs and the slowness. In [1, 28] the problem of tensor completion using the
t-product and the c-product with regularized total variation gave good
results. The problem of those methods is the fact that they are applied only
to third-order tensors. Next, we propose to extend those methods to high order
tensors using the $*_{\mathcal{L}}$-product we defined in the preceeding
section for tensors of order greater than three.
The main optimization problem that solves the problem of tensor completion
consists in finding a low-rank tensor that contains the main information (the
known pixels), which depends on the definition of the rank that we will
consider. In our proposed method, we will consider the average rank of a
tensor. Thus our main optimization problem is given as follows

$\displaystyle\underset{\mathcal{X}\in\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}\times\dots\times I_{N}}}{\arg\,\min}$ $\displaystyle\,{\tt
rank}_{a}\left(\mathcal{X}\right)$ (61) $\displaystyle{\tt s.t}$
$\displaystyle\mathcal{P}_{\Omega}\left(\mathcal{X}\right)=\mathcal{P}_{\Omega}\left(\mathcal{M}\right),$

where $\mathcal{X}$ is the underlying tensor,
$\mathcal{M}\in\mathbb{K}^{I_{1}\times I_{2}}_{I_{3}\times\dots\times I_{N}}$
is the observed tensor, $\Omega$ is the set of the known pixels and
$\mathcal{P}_{\Omega}$ is the projection operator that copy the values of the
pixels onto $\Omega$. However, the optimization problem (61) is NP-hard [23].
It is known that convex optimization problem are the easiest optimization
problems to solve and for this reason we will use the approximation of the
function average rank given in Theorem 23. Our optimization problem is
transformed to the following one

$\displaystyle\underset{\mathcal{X}\in\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}\times\dots\times I_{N}}}{\min}$
$\displaystyle\,\left\|\mathcal{X}\right\|_{*}$ (62) $\displaystyle{\tt s.t}$
$\displaystyle\mathcal{P}_{\Omega}\left(\mathcal{X}\right)=\mathcal{P}_{\Omega}\left(\mathcal{M}\right),$

where $\left\|.\right\|_{*}$ is the tensor nuclear norm in Theorem 22. The
main techniques for solving (62) is the Proximal Gradient Algorithm (PGA). The
problem (62) will be solved iteratively as

$\mathcal{X}^{k+1}=\underset{\mathcal{X}\in\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}\times\dots\times
I_{N}}}{\arg\,\min}\mu\left\|\mathcal{X}\right\|_{*}+\dfrac{1}{2}\left\|\mathcal{X}-\mathcal{G}^{k}\right\|_{F}^{2},$
(63)

where

$\displaystyle\mathcal{G}^{k}$
$\displaystyle=\mathcal{Y}^{k}-\left(\mathcal{P}_{\Omega}(\mathcal{Y}^{k})-\mathcal{P}_{\Omega}(\mathcal{M})\right)=\mathcal{P}_{\Omega^{c}}(\mathcal{Y}^{k})-\mathcal{P}_{\Omega}(\mathcal{M}),\,and$
(64) $\displaystyle\mathcal{Y}^{k}$
$\displaystyle=\mathcal{X}^{k}+\dfrac{t_{k-1}-1}{t_{k}}\left(\mathcal{X}^{k}-\mathcal{X}^{k-1}\right).$
(65)

Notice that in this case, the function $f$ involved in (10) is given by
$f(\mathcal{X})=\dfrac{1}{2}\left\|\mathcal{P}_{\Omega}(\mathcal{X})-\mathcal{P}_{\Omega}(\mathcal{M})\right\|_{F}^{2}$
with $\nabla
f(\mathcal{X})=\mathcal{P}_{\Omega}(\mathcal{X})-\mathcal{P}_{\Omega}(\mathcal{M})$.
Therefore, the Lipshitz constant $l_{f}$ of $\nabla f$ is equal to $1$. It is
clear from Theorem 25, that the solution of the problem (63) is the t-svt of
the tensor $\mathcal{G}^{k}$, i.e.,

$\mathcal{X}^{k+1}=\mathcal{D}_{\mu}\left(\mathcal{G}^{k}\right).$ (66)

THe following algorithm summarizes all the steps of the proposed method

Algorithm 5 Tensor completion using tensor nuclear norm by PGA.

1: Inputs: $\mathcal{M}\in\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}\times\dots\times I_{N}}$, $tol$, $\nu$, $\mu_{0}$.

2: Initialize: $\mathcal{X}^{0}=\mathcal{X}^{1}=0$, $t_{0}=t_{1}=0$,
$\bar{\mu}=\nu\mu_{0}$, $itermax=100$, k=1.

3: while not converged do

4: Update $\mathcal{\mathcal{Y}}^{k}$ from (65).

5: Update $\mathcal{\mathcal{G}}^{k}$ from (64).

6: Update $\mathcal{X}^{k+1}$ from (66).

7: $t_{k+1}=\dfrac{1+\sqrt{4t_{k}^{2}+1}}{2}$.

8: $\mu_{k+1}=max\left(\nu\mu_{k},\,\bar{\mu}\right)$.

9: end while

In Table 1, we give the cost of the different tensor operations (
$*_{\mathcal{L}}$-product, $*_{\mathcal{L}}$-svd and $*_{\mathcal{L}}$-svt)
when using FFT or DCT.

$*_{\mathcal{L}}$-product of $\mathcal{X}\in\mathbb{K}^{I_{1}\times
n}_{I_{3}\times\dots\times I_{N}}$ and $\mathcal{Y}\in\mathbb{K}^{n\times
I_{2}}_{I_{3}\times\dots\times I_{N}}$
---
FFT |  DCT
$O\left((PI_{1}n+PI_{2}n+PI_{1}I_{2})\sum_{i=3}^{N}\log(I_{i})\right)$ | $O\left((PI_{1}n+PI_{2}n+PI_{1}I_{2})\sum_{i=3}^{N}\log(I_{i})\right)$
$+4O\left(PI_{1}n^{2}I_{2}\right)$ | $+O\left(PI_{1}n^{2}I_{2}\right)$
$*_{\mathcal{L}}$-svd of $\mathcal{X}\in\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}\times\dots\times I_{N}}$
FFT |  DCT
$O\left((2PI_{1}I_{2}+PI_{2}^{2}+PI_{1}^{2})\sum_{i=3}^{N}\log(I_{i})\right)$ | $O\left((2PI_{1}I_{2}+PI_{2}^{2}+PI_{1}^{2})\sum_{i=3}^{N}\log(I_{i})\right)$
$+2O\left(P\min(I_{1}^{2}I_{2},I_{1}I_{2}^{2})\right)$ | $+O\left(P\min(I_{1}^{2}I_{2},I_{1}I_{2}^{2})\right)$
$*_{\mathcal{L}}$-svt of $\mathcal{X}\in\mathbb{K}^{I_{1}\times
I_{2}}_{I_{3}\times\dots\times I_{N}}$
FFT |  DCT
$O\left((2PI_{1}I_{2}+PI_{2}^{2}+PI_{1}^{2})\sum_{i=3}^{N}\log(I_{i})\right)$ | $O\left((2PI_{1}I_{2}+PI_{2}^{2}+PI_{1}^{2})\sum_{i=3}^{N}\log(I_{i})\right)$
$+2O\left(P\min(I_{1}^{2}I_{2},I_{1}I_{2}^{2})\right)+8O\left(PI_{1}n^{2}I_{2}\right)$ | $+O\left(P\min(I_{1}^{2}I_{2},I_{1}I_{2}^{2})\right)+2O\left(PI_{1}n^{2}I_{2}\right)$

Table 1: The cost of computing the $*_{\mathcal{L}}$-product,
$*_{\mathcal{L}}$-svd and the $*_{\mathcal{L}}$-svt by using Fourier and
cosine transforms.

##  5 Numerical experiments

In this section we test the performance of our algorithms for high order
tensor completion using the Fourier and Cosine for thethe operator
$\mathcal{L}$ and we will compare the obtained results with those ontained by
some existing known algorithms on color videos. In Subsection 5.1, the tests
were performed with Matlab 2018a, on an Intel i5 laptop with 16 Go of memory,
and in Subsection 5.2 we use codes with Python on a machine that uses a CPU of
type Intel Xeon Gold 6152 with a frequency from 2.1Ghz to 3.7 GHz and a GPU of
type NVIDIA Tesla Pascal 40. All the tests are computed using a single core.
The quality of the obtained data can be computed by the peak signal-to-noise-
ration (PSNR) defined by

$PSNR=10\,\log_{10}\dfrac{Max_{\mathcal{A}_{obt}}^{2}}{\left\|\mathcal{A}_{obt}-\mathcal{A}_{ori}\right\|_{F}^{2}},$
(67)

and the relative squared error (RSE) given by

$RSE=\dfrac{\left\|\mathcal{A}_{ori}-\mathcal{A}_{obt}\right\|_{F}^{2}}{\left\|\mathcal{A}_{obt}\right\|_{F}^{2}},$
(68)

where $\mathcal{A}_{ori}$ is the original tensor, $\mathcal{A}_{obt}$ is the
obtained recovered tensor and $Max_{\mathcal{A}_{obt}}$ is the maximum pixel
of the recovered tensor. The quality of the recovered data is good when the
value of RSE is small and the value of PSNR is high. In our experiments of
tensor completion, we use $M_{i}=F_{i}$ for Fourier transform and
$M_{i}=C_{i}$ for Cosine transform, where the transformed data
$\mathcal{L}(\mathcal{X})$ for an $N^{th}$-order tensor is given by ${\tt
fft}(\mathcal{X},[\,],3)$ and ${\tt dct}(\mathcal{X},[\,],3)$, respectivly.
The parameters $\mu_{0}$ and $\nu$ stated the algorithm of completion have to
be fixed. We set $\nu=0.9$ and $\mu_{0}=\nu\left\|\mathcal{M}\right\|_{F}$,
where the stopping criterion convergence of this algorithm is as follows

$\dfrac{\left\|\mathcal{Y}^{k}-\mathcal{X}^{k+1}\right\|_{F}}{\left\|\mathcal{X}^{k+1}\right\|_{F}}\leq
10^{-4}.$ (69)

Figure 1 shows the data tests used in our experiments. The video of xylophone
is available from Matlab, and the videos Akiyo and News are available from
111urlhttp://trace.eas.asu.edu/yuv/.

|  |
---|---|---
xylophone | Akiyo | news
Figure 1: The $20^{th}$ bound of the color videos of the test data.

In Table 2, we give the size of the different used color videos tests.

name | size
---|---
xylophone | $240\times 320\times 3\times 141$
Akiyo | $144\times 176\times 3\times 300$
news | $144\times 176\times 3\times 300$

Table 2: The name of the test data (color videos) and their sizes.

###  5.1 Tensor completion for color videos

In this part we show the obtained results of our algorithms TNN-PGA-F and TNN-
PGA-C on fourth-order tensors (color videos) and compare them with other ones
such as the methods named Tmac [22] and HaLRTC [26]. The comparison will be in
terms of the efficiency and executing times by comparing the values of RSE and
PSNR, the number of iterations the required cpu-time. In Figure 2 we show the
$20^{th}$ bound of the color videos xylophone, Akiyo and news with only $5\%$
of the original data. In Figure 3, we show the $20^{th}$ bound of the
recovered data for each color video obtained by by the algorithms HaLRTC,
Tmac, TNN-PGA-F and TNN-PGA-C.

|  |
---|---|---
xylophone | akiyo | news
Figure 2: The $20^{th}$ bound of color videos of the test data with only $5\%$
of the original data, i.e., $sr=0.05$.

|  |  |
---|---|---|---

|  |  |
---|---|---|---

|  |  |
---|---|---|---
HaLRTC | TMac | TNN-PGA-F | TNN-PGA-C

Figure 3: The $20^{th}$ bound of the recovered color videos obtained from the
algorithms HaLRTC, Tmac, TNN-PGA-F and TNN-PGA-C for $sr=0.05$.

In Table 3 we report the values of RSE, PSNR, the number of iterations and the
CPU time required by the algorithms HaLRTC, Tmac,TNN-PGA-F and TNN-PGA-C.
Figure 4 shows the curves representing the evolution of the RSE and PSNR
versus the iteration number for Akiyo-video with $sr=0.05$. In Figure 5 we
give the values of the first $1000$ pixels of the recovered data of Akiyo-
video for $sr=0.05$. Figure 6, shows the values of RSE and PSNR for each bound
of the xylophone-video for $sr=0.1$.

sr | $0.05$ | $0.1$
---|---|---
Video | Algorithm | RSE | PSNR | Iteration | time | RSE | PSNR | Iteration | time
Akiyo | HaLRTC | 0.3741 | 15.49 | 163 | 2464.8 | 0.2409 | 19.31 | 94 | 1632.4
Tmca | 0.0879 | 30.25 | 343 | 991.0 | 0.0824 | 30.91 | 210 | 587.6
TNN-PGA-F | 0.0497 | 33.84 | 73 | 619.8 | 0.0356 | 36.65 | 70 | 640.8
TNN-PGA-C | 0.0496 | 33.91 | 99 | 571.2 | 0.0350 | 36.80 | 76 | 475.1
xylophone | HaLRTC | 0.2732 | 17.49 | 163 | 3247.3 | 0.1910 | 20.60 | 96 | 2012.0
Tmca | 0.1148 | 26.96 | 376 | 3354.7 | 0.1066 | 27.43 | 347 | 1703.7
TNN-PGA-F | 0.0888 | 28.87 | 73 | 1054.2 | 0.0680 | 30.77 | 73 | 1041.4
TNN-PGA-C | 0.0889 | 28.38 | 90 | 906.0 | 0.0660 | 30.90 | 79 | 760.9
news | HaLRTC | 0.4951 | 14.46 | 142 | 1876.3 | 0.3694 | 17.01 | 89 | 1038.2
Tmca | 0.1344 | 27.95 | 1157 | 4480.0 | 0.1171 | 29.36 | 415 | 1752.6
TNN-PGA-F | 0.1058 | 29.07 | 71 | 621.9 | 0.0796 | 31.53 | 69 | 573.4
TNN-PGA-C | 0.1156 | 28.20 | 99 | 614.0 | 0.0838 | 31.15 | 77 | 475.5

Table 3: The values of RSE, PSNR, the number of iteration and the time
required by the algorithms HaLRTC, Tmac, TNN-PGA-F and TNN-PGA-C for $sr=0.05$
and $sr=0.1$.

|
---|---
Figure 4: The evolution of the RSE and the PSNR values on each iteration of
the algorithms HaLRTC, Tmac, TNN-PGA-F and TNN-PGA-C for the video of Akiyo
with $sr=0.05$.
---

Figure 5: Comparison of the first $1000^{th}$ pixels of the video Akiyo obtained by the Algorithms TNN-PGA-F and TNN-PGA-F with $sr=0.05$ by the original data of the video. |
---|---

Figure 6: The values of RSE and PSNR ones each bound of the recovered video
xylophone for $sr=0.1$.

From Figures 3 we can see that our algorithms return very good results. Table
3 confirms this fact showing an advantage for the cosine transform as compared
to the Fourier transform. Figure 4 shows that the curves obtained by our
algorithms decreases (for RSE) and increases (for PSNR) quickly towards the
minimum and the maximum value, respectively. Figure 5 shows that the values of
the recovred data obtained by the proposed two algorithms are very close to
the original data. Figure 6 shows the efficiency of our algorithms as compared
to other ones for each bound of the used video.

###  5.2 Porting the python code for information completion to GPU using CuPy

Parallel computation can be very important in high performance computing due
to the limit of the use of a single core. As a consequence, the latest CPU
manufacturers compete to have the most cores on a single CPU. For highly
parallel problems that can benefit from more cores using a GPU (Graphics
Processing Unit), it is the best approach as GPU sacrifices memory for more
cores per unit and this leads to having a massively parallel capabilities over
CPU. The mainly use of a GPU is to rapidly manipulate and alter memory to
accelerate the creation of images.
In our tests, all the GPU accelerated libraries utilize CUDA toolkit libraries
which is a parallel computing platform and an API that allows interaction with
a GPU in order to perform general purpose processing. That goes beyond just
image data manipulation allowing a more general approach called GP-GPU which
stands for general purpose computing on graphics processing units. We used
Cupy to accelerate some part of the CPU code. CuPy is an open-source library
with NumPy syntax that increases speed by doing matrix operations on NVIDIA
GPUs. It is accelerated with the CUDA platform from NVIDIA and also uses CUDA-
related libraries, including cuBLAS, cuDNN, cuRAND, cuSOLVER,cuSPARSE, and
NCCL, to make full use of the GPU architecture. CuPy’s interface is highly
compatible with NumPy and in most cases it can be used as a drop-in
replacement that can easily integrated in already existing CPU code to boost
the performance without much code changes.

####  5.2.1 Porting the python code to GPU using CuPy

Cupy provides an easy way to port a Python code using Numpy and Scipy by
accelerating them using GPU. The porting process can be simple by replacing
some Numpy/Scipy functions by their equivalent in Cupy.

Code of tsvt by fft and dct:

  * •

Before using CuPy

  * •

After using Cupy

For the case of of tsvt-cosine by CPU we will use the same steps by changing
fftn and ifftn by dctn and idctn, respectively This code provides a good
example on how we can introduce Cupy to an already existing peace of code by
identifying the heavy work functions and see if they have a Cupy equivalent.
The most time consuming functions are svd, fftn and ifftn. We can see in the
code below that cupy provides a GPU-accelerated implementation of those
functions . For the algorithm that uses the DCT function instead of the FFT,
the problem was more difficult since Cupy don’t support the parallelizable
version.

####  5.2.2 Numerical experiments of the problem of completion

In Table 4 we give the size of all the data test used in our experiments
(color videos: fourth-order tensors).

Name | size
---|---
xylophone | $240\times 320\times 3\times 30$
car | $1920\times 1080\times 3\times 30$
Mgrass | $2160\times 4096\times 3\times 30$
notes | $3840\times 2160\times 3\times 30$

Table 4: The size of all the data used in the experiments.

In this part we give the results of our codes of completion (TNN-PGA-F and
TNN-PGA-C) by using CPU and GPU computation. In the next we denote by PGA-F
and PGA-C the codes using CPU and by PGA-F-GPU and PGA-C-GPU those using GPU.
In Figures 7 and 8, we compare the evolution of the RSE and the error during
the execution of the codes by CPU and GPU for two different videos ’mglass’
and ’notes’ with two values of $sr$, $sr=0.05$ and $sr=0.1$. In Figure 9 we
give an histogramme representing the required time of PGA-F-CPU, PGA-F-GPU,
PGA-C-CPU and PGA-C-GPU for $sr=0.1$.

|
---|---
Figure 7: The evolution of RSE and the error on each iteration for mglass starting with $5\%$ of the original data |
---|---
Figure 8: The evolution of RSE and the error on each iteration for notes
starting with $10\%$ of the original data

Figures 7 and 8 show that the RSE and the error does not change when using the
Cupy function as they are almost identical to the Numpy and Scipy ones.

---
Figure 9: time comparison for data starting with $10\%$ of the original data.

In Table 5 we show the speed-up between the sequential and parallel
computations. We reported the average $\tau$ defined by using the following
formula

$\tau=\dfrac{cpu\\_time}{gpu\\_time}.$ (70)

| xylohone | car | mglass | notes
---|---|---|---|---
PGA-C | $0.77$ | $1.68$ | $1.5$ | $1.80$
PGA-F | $0.86$ | $2.11$ | $2.47$ | $2.89$

Table 5: Speed up percentage for the data set that started from $10\%$ from
the original data.

As shown in Table 5, the obtained speed-up values show how much we can boost
the performance of our codes by using the GPU. When the data set is small, as
in the case for xylophone, there is no need to use GPU becuase in that case
the returned cpu-time is smaller than the one optained by GPU. This
performance can be explained by the fact that for small problems, the transfer
of the data in parallel computation, requires a significant time compared the
classical computation for which no need of transfert data is needed. When the
data becomes larger and the targeted accelerated function is taking a
significant time from the total runtime we can see a big speed up using the
GPU accelerated functions as kernels.

##  6 Conclusion

In this paper we presented a new tensor-tensor product for high orders. Using
this product, we defined a new high -order SVD and some related properties. We
gave some theoretical results for the tensor product. We used this tensor
product for tensor completion using the proximal gradient algorithm. In the
numerical section, we showed some test on color-videos and used GPU
computation to get fast computation . The presented numerical experiments show
the efficiency of our proposed algorithms.

## References

  * [1] S. Aeron, G. Ely, N. Hoa, M. Kilmer, Z. Zhang, Novel methods for multilinear data completion and de-noising based on tensor-SVD. J. Proceedings of the IEEE conference on computer vision and pattern recognition, 3842–3849 (2014).
  * [2] S. Aeron, E. Kernfeld, M. Kilmer, Tensor–tensor products with invertible linear transforms, Linear Algebra and its Applications, 485, 545–570 (2015).
  * [3] H. Avron, L. Horesh, M. Kilmer, E. Newman, Tensor-tensor algebra for optimal representation and compression of multiway data, Proceedings of the National Academy of Sciences, 118, 28 (2021).
  * [4] B. W. Barder, T. Kolda, Tensor decompositions and applications, SIAM review, 51 (3), 455–500 (2009).
  * [5] A. Beck, M. Teboulle, A fast iterative shrinkage-thresholding algorithm for linear inverse problems, SIAM journal on imaging sciences, 2 (1), 183–202 (2009).
  * [6] S. P. Boyd, M. Fazel, H. Hindi, A rank minimization heuristic with application to minimum order system approximation. Proceedings of the 2001 American Control Conference.(Cat. No. 01CH37148). 6, 4734–4739 (2001).
  * [7] K. Braman, N. Hao, R. C. Hoover, M. E. Kilmer, Third-order tensors as operators on matrices: A theoretical and computational framework with applications in imaging, SIAM Journal on Matrix Analysis and Applications, 34, (1) 148–172 (2013).
  * [8] J. A. Bengua, M. N. Do, H. N. Phien, H. D. Tuan, Efficient tensor completion for color image and video recovery: Low-rank tensor train, IEEE Transactions on Image Processing, 26, (5) 2466–2479 (2017).
  * [9] A. H. Bentbib, A. El Hachimi, K. Jbilou, A. Ratnani, A Tensor Regularized Nuclear Norm Method for Image and Video Completion, Journal of Optimization Theory and Applications, 1–25 (2021).
  * [10] A. H. Bentbib, A. El Hachimi, K. Jbilou, A. Ratnani, Fast multidimensional completion and principal component analysis methods via the cosine product.
  * [11] J. F. Cai, E. J. Candès, Z. Shen, A singular value thresholding algorithm for matrix completion, SIAM Journal on optimization, 20 (4), 1956–1982 (2010).
  * [12] C. D. Carla, M. E. Kilmer, Factorization strategies for third-order tensors. Linear Algebra and its Applications, 435, 641–658 (2011).
  * [13] Y. Chen, J. Feng, H.Lin, W. Liu, C. Lu, S. Yan, Tensor robust principal component analysis with a new tensor nuclear norm. IEEE transactions on pattern analysis and machine intelligence, 42, 925–938 (2019).
  * [14] M. Chen, A. Ganesh, Z. Lin, Y. Ma, L. Wu, J. Wright, Fast convex optimization algorithms for exact recovery of a corrupted low-rank matrix,Coordinated Science Laboratory Report no. UILU-ENG-09-2214, DC-246, (2009).
  * [15] G.H. Golub, C.F. Van Loan, Matrix Computations, 4th edition, Johns Hopkins University Press, Baltimore, MD, 2013.
  * [16] L. De Lathauwer, Decompositions of a higher-order tensor in block terms — Part II: Definitions and uniqueness, SIAM J. Matrix Anal. Appl., 30, 1033–1066, (2008).
  * [17] L. De Lathauwer, B. De Moor, J. Vandewalle, A multilinear singular value decomposition, SIAM journal on Matrix Analysis and Applications, (21) 4, 1253–1278 (2000).
  * [18] M. Ding, T. -Z. Huang, T. -Y. Ji, J. -H. Yang, X. -L. Zhao, Low-rank tensor completion using matrix factorization based on tensor train rank and total variation. Journal of Scientific Computing, 81, 941–964 (2019).
  * [19] A.El Ichi, K. Jbilou and R. Sadaka, On tensor tubal-Krylov subspace methods, Linear and Multilinear Algebra, 2021, DOI: 10.1080/03081087.2021.1999381
  * [20] M. El Guide , A. El Ichi, K. Jbilou, R. Sadaka,: On tensor GMRES and Golub-Kahan methods via the T-product for color image processing. The Electronic Journal of Linear Algebra. 37, 524–543 (2021).
  * [21] D. Goldfarb, Z. Qin, Robust low-rank tensor recovery: Models and algorithms, SIAM Journal onMatrix Analysis and Applications, 35 (1), 225–253 (2014).
  * [22] R. Hao, Z. Su, Y. Xu, W. Yin, Parallel matrix factorization for low-rank tensor completion, arXiv preprint arXiv:1312.1254, (2013).
  * [23] C. J. Hillar, L. -H. Lim, Most tensor problems are NP-hard, Journal of the ACM (JACM), 60, 1–39 (2013).
  * [24] T. -Y. Ji, T. -Z. Huang, G. Liu, T. -H. Ma, X. -L. Zhao, Tensor completion using total variation and low-rank matrix factorization, Information Sciences, 326, 243–257 (2016).
  * [25] M.E. Kilmer, K. Braman, N. Hao, R.C. Hoover, Third-order tensors as operators on matrices: a theoretical and computational framework with applications in imaging, SIAM J. Matrix Anal. Appl. 34 (2013) 148–172.
  * [26] J. Liu, P. Musialski, P. Wonka, J. Ye, Tensor completion for estimating missing values in visual data, IEEE transactions on pattern analysis and machine intelligence, 35 (1), 208–220 (2012).
  * [27] C. Martin, B. LaRue, R. Shafer, An order-p tensor factorization with applications in imaging, SIAM Journal on Scientific Computing, 35, 1 (2013).
  * [28] M. Ng ,W. -H. Xu, X. -L. Zhao, A fast algorithm for cosine transform based tensor singular valuedecomposition, arXiv preprint arXiv:1902.03070, (2019).