Post-Processing of High-Dimensional Data

by   Hermann G. Matthies, et al.
RWTH Aachen University

Scientific computations or measurements may result in huge volumes of data. Often these can be thought of representing a real-valued function on a high-dimensional domain, and can be conceptually arranged in the format of a tensor of high degree in some truncated or lossy compressed format. We look at some common post-processing tasks which are not obvious in the compressed format, as such huge data sets can not be stored in their entirety, and the value of an element is not readily accessible through simple look-up. The tasks we consider are finding the location of maximum or minimum, or minimum and maximum of a function of the data, or finding the indices of all elements in some interval --- i.e. level sets, the number of elements with a value in such a level set, the probability of an element being in a particular level set, and the mean and variance of the total collection. The algorithms to be described are fixed point iterations of particular functions of the tensor, which will then exhibit the desired result. For this, the data is considered as an element of a high degree tensor space, although in an abstract sense, the algorithms are independent of the representation of the data as a tensor. All that we require is that the data can be considered as an element of an associative, commutative algebra with an inner product. Such an algebra is isomorphic to a commutative sub-algebra of the usual matrix algebra, allowing the use of matrix algorithms to accomplish the mentioned tasks. We allow the actual computational representation to be a lossy compression, and we allow the algebra operations to be performed in an approximate fashion, so as to maintain a high compression level. One such example which we address explicitly is the representation of data as a tensor with compression in the form of a low-rank representation.



There are no comments yet.



Computing f-Divergences and Distances of High-Dimensional Probability Density Functions – Low-Rank Tensor Approximations

Very often, in the course of uncertainty quantification tasks or data an...

Multi-resolution Low-rank Tensor Formats

We describe a simple, black-box compression format for tensors with a mu...

Tensor-Tensor Products for Optimal Representation and Compression

In this era of big data, data analytics and machine learning, it is impe...

Tensor-train approximation of the chemical master equation and its application for parameter inference

In this work, we perform Bayesian inference tasks for the chemical maste...

Tensor-based EDMD for the Koopman analysis of high-dimensional systems

Recent years have seen rapid advances in the data-driven analysis of dyn...

Improving Matrix-vector Multiplication via Lossless Grammar-Compressed Matrices

As nowadays Machine Learning (ML) techniques are generating huge data co...

ALTO: Adaptive Linearized Storage of Sparse Tensors

The analysis of high-dimensional sparse data is becoming increasingly po...
This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

1 Introduction

Many scientific and engineering computations or measurements, as well as economic or financial applications, produce large volumes of data. For simplicity, assume that we have computed or observed large quantities of one real-valued variate. Assume further that the amount of data is so large that it can not be held or stored in its entirety and has to be compressed in some way.

Frequently, the task arises for example to find the location with the maximum value of the data set, or to find the locations of all values which lie in a given interval — we call this a level set. Other, similar tasks we consider are finding the number of values in a given interval, the probability of being in a given interval, or finding the mean or the variance of the data set. These tasks are trivial for small data sets, typically performed by inspecting each datum. But if we consider truly huge amounts of data which can not be stored in full because this would exceed the capacity of any storage device, but only in a compressed manner, and where additionally it may not be possible due to time constraints to inspect each element, these tasks are not trivial any more. This is due to the fact that in the compressed representation the data values are normally not directly accessible and require additional processing. Such compression will in general be “lossy”, so that not each value can be restored exactly, but only up to a certain reconstruction error.

Assuming the data as an element of some set , the algorithms are independent of the representation of the data, as well as from the compression and reconstruction technique used, subject only to the possibility of approximately performing the operations of an Euclidean associative commutative algebra. This means that we assume that

is a vector space with an inner product, and additionally an associative commutative bilinear multiplication, making it into an algebra. The simplest example of such a structure is to envisage the data points as a vector

with an appropriate , reflecting the amount of data. The vector space operations are clear, as is the canonical Euclidean product, and for the commutative multiplication consider the point-wise or Hadamard product, i.e. the component-wise multiplication of two vectors. We shall allow that all these operations are performed only approximately, in order to maintain a high compression level.

Large volumes of data, especially when they can be thought of as samples or discrete values of some real-valued function on a high-dimensional space, can often be arranged in form of a tensor [40, 41, 8, 31]

. This offers the possibility to use approximate or compressed tensor representations. Here we will especially show one possibility, namely low-rank representations, which generalise the truncated singular value decomposition for matrices. For the sake of completeness, we shall show the possible implementation of the above mentioned algebraic operations for some of the more common low-rank formats, and the approximated, compressed, or truncated representation of their results; see also


The proposed algorithms are iterative in nature, and the convergence tolerance can be adapted to the reconstruction error. The basic idea for the algorithms, which operate only on the algebraic structure, is an iterative scheme which converges to a result which solves the desired problem in some way. Such iterations typically destroy the compressed representation, so they have to be combined with re-compression or truncation.

1.1 An example

Let us give an example which motivates much of the following formulation and development. Assume that we are interested in the time evolution of some system, described by


where is in some Hilbert space and is some parameter dependent operator; in particular could be some parameter-dependent differential operator, for example


where is a domain, is the time window of interest, is a parameterised random tensor field dependent on a random parameter in some probability space with probability measure , and one may take [43, 46].

For each one may thus seek for solutions in , where we have set . On the other hand, assume that for fixed

we are looking at the random variable

which we for simplicity assume to have finite variance, i.e. . Taking into account the parametric dependence, we are hence looking for a function which is defined on ; we are thus looking for a solution in . This applies equally well to the abstract eq:XXV, as on the right hand side the operator depends on , so will the solution to eq:XXV, and it will lie in a similar tensor product.

Observe further that if the probability measure on in eq:XXV or eq:XXVI is a product measure on , the space can be split further , with . Thus in total the solution to both eq:XXV and eq:XXVI may be viewed as an element of a high order or degree tensor space [43, 46, 48]


If we think of samples of at space points , time instances and parameter values , with a multi-index with , the samples of the solution


form a high order or degree tensor [31].

In many chemical applications one has , and is few hundreds or even thousands [6]. Some other high-dimensional problems from chemistry and physics like Hartree-Fock-, Schrödinger-, or Master-equations, and low-rank tensor methods of their solution are for example considered in [32, 37, 39, 11], see also the references therein. Another example of large volumes of high-dimensional data are satellite data. Satellites collect data over a very large areas (e.g. the data collected by the National Center for Atmospheric Research (USA) [29]. Big data can also come from a computer simulator codes such as a solution of a multi-parametric equation, e.g. weather research and forecasting, and climate models [23]. Also OilGas companies daily collect sensor data from multiple sources. Another source for huge volumes of data are high-energy particle accelerators like CERN [5].

For the sake of simplicity, we will treat all arguments of the function in eq:XXaX equally, denoting them in that case by . So just consider a real-valued function and its samples


which form a tensor of order or degree representing a total of values, and where for the sake of convenience multi-indices have been introduced.

As one may see, a tensor can be simply defined as a high-order matrix or multi-index array, where multi-indices are used instead of indices. As an example, assume that we have data points. Now a vector is just a tensor of first order, but one may reshape the data, say into a matrix — a tensor of 2nd order. Reshaping further, one may take — a tensor of 4th order — as well as other combinations, such as — an array of size (16 times) — which is a tensor of 16th order; and finally all the way to — a tensor of order 32.

The tensors obtained in this way contain not only rows and columns, but also slices and fibres [40, 41, 8, 31]. These slices and fibres can be analysed for linear dependencies, super symmetry, or sparsity, and may result in a strong data compression. To have a first glimpse of possible compression techniques, assume that in the above example the data has been stored in the matrix , where , and consider its singular value decomposition (SVD) , where is the diagonal matrix of singular values , assumed arranged by decreasing value, and , collect the right and left singular vectors. Often there is a number , such that for some small one has for all . Then one can formulate a compressed or truncated version

where , , , and , . This reduced SVD is a special case of eq:CPI. For the sake of a concrete example, assume that . One may observe that the required storage has been reduced from to for the vectors .

The set of possible objects we want to work with will be denoted by . We assume that this set is equipped with the operations of an associative and commutative real algebra, and carries an inner product. It is a well known result [60] that such algebras — modulo some technical details — are isomorphic via the Gel’fand representation to a function algebra — to be more precise continuous real valued functions on a compact set. Under point-wise addition and multiplication by scalars, such functions clearly form a vector space, and if one includes point-wise multiplication of functions, they form an associative and commutative real algebra. In our case this is simply the function algebra , which is obviously the same as . The advantage of this abstract algebraic formulation is not only that it applies to any kind of data representation on which one may define these algebraic operations, but that one may use algorithms [35] which have been developed for the algebra of real matrices .

In our concrete examples we work directly with the data represented as tensors . Even without the identification with , it is obvious that this is naturally a vector space. As also , it is clear that can be equipped with the canonical Euclidean inner product from . The associative and commutative algebra product comes from the identification with the function algebra , i.e. the point-wise product, which is also known as the Hadamard product [31]. Hence it is clear that for data it is not difficult at all to define the algebraic operations, but rather how to perform them when the data is in a compressed format.

To simplify later notation on the number of operations and amount of storage needed, we will often make the assumption that , so that the tensor in eq:XXbX represents values. As already mentioned, our example of data compression is based on low-rank approximations to elements in . Although low-rank tensor data formats and techniques are almost unavoidable if working with large high-dimensional data sets, we would like to stress that the algorithms presented here are formulated purely in terms of the abstract algebraic structure, and are thus independent of the particular representation.

Whereas a general element hence has terms, a compressed representation — some low-rank versions of which will be discussed in more detail in tensor-rep — will have significantly fewer terms. For example, the CP-decomposition (canonical polyadic) representation, truncated to terms,


If the rank is reasonably small compared to independent of and , then we have an approximation with much less storage, which also only depends on the product of and and is not exponential in the dimension . But now the maximum can not be easily looked up; to get a particular element, the expression eq:CPI has to be evaluated for that index combination. In the following we shall assume that we work with approximations such as in eq:CPI which need much less storage. One has to make sure then that the algebraic operations of the Hadamard algebra structure on can be performed efficiently in an approximative manner, so that they have much lower complexity than the obvious , although the single elements of which are usually needed for point-wise operations — and post-processing, see post-proc-intro — are not directly available. Thus the motivating factors for applying compression, and in particular low-rank tensor techniques, include the following:

  • The storage cost is reduced, depending on the tensor format, from to , or to , where .

  • The low-rank tensor approximation is relatively new, but already a well-studied technique with free software libraries available. Various low-rank formats are available.

  • The approximation accuracy is fully controlled by the tensor rank. The full rank gives an exact representation.

  • Even more complicated operations like the Fourier transform can be performed efficiently in low-rank format. The basic fast Fourier transform on

    would have complexity . These low-rank techniques can either be combined with the fast Fourier transform giving a complexity of , or can even be further accelerated at the price of an additional approximation error yielding a superfast Fourier transform [49, 10].

On the other hand, general limitations of a compression technique are that

  • it could be time consuming to compute a compression, or in particular a low-rank tensor decomposition;

  • with sampling, it requires an axes-parallel mesh;

  • although many functions have a low-rank representation, in practice only some theoretical estimates exist.

But the fact still remains that there are situations where storage of all items is not feasible, and some kind of compression has to be employed.

1.2 Post-processing tasks

This work is about exploiting the compression together with the structure of a commutative algebra and the ability to perform the algebraic operations in the compressed format, at least approximately. The tensor product structure which appears in eq:XXIX is something which allows efficient calculations to be performed on a sample of the solution like eq:XXaX or eq:XXbX with the help of expressions such as eq:CPI or other compression techniques.

This tensor product structure—in this case multiple tensor product structure—is typical for such parametric problems [47]. What is often desired, is a representation which allows for the approximate evaluation of the state of eq:XXV or eq:XXVI as function of without actually solving the system again. Furthermore, one would like this representation to be inexpensive to evaluate, and for it to be convenient for certain post-processing tasks, for example like finding the minimum or maximum value over some or all parameter values.

The tasks we consider here are

  • finding the location and value of maxima or minima,

  • finding the location and value of maxima or minima of a function of ,

  • finding the value in the tensor closest to some given number,

  • finding level sets, i.e. the indices where the value is between two levels,

  • finding the number of indices between two levels,

  • computing the probability of being between two levels,

  • computing the mean and the variance.

As an example, consider finding the maximum value. A naïve approach to compute the maximum would be to visit and inspect each element, but then the number of visits is , exponential in , which may be unacceptable for high

. Such phenomena when the total complexity/storage cost depends exponentially on the problem dimension have been called the curse of dimensionality


The idea for such iterative post-processing algorithms, in the form of finding the maximum, was first presented in [14] for low-rank tensor approximations. Some further post-processing tasks for such low-rank tensor representations were investigated in [20]. To give an example from [14] of the possible savings in space and time, assume that one looks at the solution of an example of eq:XXV in xmpl-motiv, or more specifically eq:XXVI in the stationary case with , a -dimensional Poisson equation:

and right-hand-side

Assume further that this is solved numerically by a standard finite-difference method with grid-points in each direction, so that the solution vector has data points in , but can also naturally be viewed as a tensor of degree in . So for full storage, one needs storage locations, and if one is looking for the maximum by inspecting each element, one would have to inspect elements.

We mention in passing that of course actually solving the discrete equation is a challenge for , but that is a different story, and for this we refer to [14] and the references therein.

# loc’s.: years [a] actual [14] time [s]
inspect. of Algorithm 2
Age of the universe: years. Number of hadrons (elementary particles) in the universe .
Table 1: Computing times (4th column) on 2 GHz dual-core CPU to find maximum.

Now assume for the sake of simplicity that one can inspect elements per second — on an ideal 2 GHz CPU with one inspection per cycle. Then for the times needed to find the maximum for the full data-set are shown in tab:find_max in the third column—assuming that it were somehow possible to store all the values indicated in the second column—whereas the actual computation with a compressed format is shown in the last and fourth column.

It is obvious that for growing for the full representation the computational complexity and the storage requirements quickly become not only unacceptable, but totally impossible to satisfy. The second and third column behave like and grow exponentially with , whereas for a low-rank representation of not only can the data be stored on a modest laptop, but for the simple Algorithm 2 — to be explained later in auxil — the computing times in the fourth column are in terms of seconds and behave like .

1.3 State of the art

This is not a discussion of the state of the art regarding general tensor formats and their low-rank approximations, in quantum physics also known as tensor networks. For the physical motivations and numerical developments from there see [64, 59, 22, 50, 4, 3]. For the mathematical and numerical view we refer to the review [41], the monographs [31, 37, 39], and to the literature survey on low-rank approximations [26]. In the following, we concentrate on the question of post-processing such data.

The idea of finding the largest element and the corresponding location of a tensor by solving an eigenvalue problem, where the matrix and the vectors are tensors given in the CP tensor format, was introduced in [14]. Additionally, the first numerical schemes to compute the point-wise inverse and the sign function were introduced, as well as the rank-truncation procedure for tensors given in the CP format. Later, in [20], these ideas were extended and applied to tensors which were obtained after discretisation and solution of elliptic PDEs with uncertain coefficients. Another group of authors, in [9], by combining the advances of the density matrix renormalisation group and the variational numerical renormalisation group methods, approximated several low-lying eigenpairs of large Hermitian matrices simultaneously in the block version of the TT format via the alternating minimisation of the block Rayleigh quotient sequentially for all TT cores.

In [13, 12], the authors suggested methods to compute the mean, the variance, and sensitivity indices in the TT train tensor format with applications to stochastic PDEs, whereas the diagonalisation of large Toeplitz or circulant matrices via combination of the fast Fourier and the CP tensor format was shown in [49].

An investigation of approximations to eigenfunctions of a certain class of elliptic operators in

by finite sums of products of functions with separated variables is the topic of [33]. Various tensor formats were used for a new class of rank-truncated iterative eigensolvers. The authors were able to reduce the computational cost from to . They demonstrated the introduced algorithms by solving large-scale spectral problems from quantum chemistry: the Schrödinger, the Hartree–Fock, and the Kohn–Sham equations in electronic structure calculations.

1.4 Outline of the paper

In the following algos the necessary material for abstract algebras is quickly reviewed. Then the algorithms and functions on the algebra used to compute the post-processing tasks outlined in post-proc-intro are formulated in an abstract fashion, independent of the representation chosen for the data. Even as the formulation uses only the abstract algebra operations, we point out and motivate what this means in terms of the Hadamard algebra. But also for the Hadamard algebra, the presentation is independent of the tensor format to be chosen. Furthermore, for the iterative algorithm — a fixed point iteration — it is discussed how the possible truncation operation influences the convergence.

In tensor-rep we present some concrete examples of compression of high-dimensional data , once it is identified with a tensor with . For such tensors of degree we discuss the canonical polyadic (CP) representation, the Tucker representation, and the tensor-train (TT) representation, and the compression in terms of low-rank approximations based on the different tensor formats. In particular, we show how the algebra operations can be carried out in the different tensor formats, the numerical effort involved, and the effect these algebraic operations have on the compression or low-rank representation. As the compression level may deteriorate, i.e. the rank of the approximation may grow, it is important that one is able to re-compress. Pointers to such methods in the literature are included as well.

The num-exp contains a few numerical examples and a discussion of their results used to illustrate the algorithms of algos and the representations from

tensor-rep. These examples are solutions of elliptic high-dimensional or parametric resp. stochastic partial differential equations, a domain where the data is naturally in a tensor format. But as already pointed out, for any data it is just a — often only conceptual —

reshape to consider it as an element of a tensor space. The conclusion is contained in concl.

2 Algorithms for post-processing

Given a data-set in some compressed tensor format, where looking up each element is not trivial, and where additionally it may not be computationally feasible to look at each element, we still want to be able to perform tasks which essentially require looking at each datum. The tasks outlined in post-proc-intro we want to consider are finding

  1. the indices or values of maxima, minima of ,

  2. the indices or values of maxima, minima of some function of ,

  3. the index or value of the datum of closest to a given number,

  4. level sets, i.e. (the indices of) all ,

  5. the number of indices such that ,

  6. the mean (sum of all items) and variance (sum of all items squared),

  7. some auxiliary functions, such as the algebraic inverse of — in our case the Hadamard inverse — and others like the function which will be seen to be needed for the above tasks in post-proc-tasks.

The above tasks and auxiliary functions will be computed by finding an iteration (mapping) for each post-processing task , such that the fixed-point of — i.e.  — is either the sought solution for the task, or computes one of the auxiliary functions.

2.1 Iteration with Truncation

When performing computations using the operations of the algebra — like the Hadamard algebraic operations on tensors in some compressed format — after the operation the compression will typically be sub-optimal, which means that the result has to be compressed again. Thus, when performing the algebraic operations for a fixed-point iteration, the compression would either get worse and worse in each iteration, or one has to re-compress or truncate the result again to a good compression level. Therefore a re-compression after each or after a number of algebraic operations may be necessary, and we thus allow that the algebraic operations are possibly only executed approximately.

In our example case the low-rank representation of tensors explained in tensor-rep acts as compression, and the rank may increase through the algebraic operations, and hence we will use truncation to low rank with error [31, 2] to avoid the problem of a deteriorating compression level.

In other words, with a general compressed representation of the computation will be a truncated or perturbed iteration [34, 48]. If we denote the general compression mapping by — meaning compression with an accuracy — the iteration map is changed from to .

The general structure of the iterative algorithms for a post-processing task or for an auxiliary function is shown in Algorithm 1. Here we want to collect some results for this kind of iteration.

1:Start with some initial compressed guess depending on task .
3:while no convergence do
4:     ;
5: the iterator may deteriorate the compression level
6:     if compression level of is too bad then
7:          ;
8: use truncation to compress with error
9:     else
10:          ;
11:     end if
13:end while
Algorithm 1 Iteration with truncation

For such a truncated or perturbed iteration as in Algorithm 1, it is known that

  1. if the iteration by is super-linearly convergent, the truncated iteration will still converge super-linearly, but finally stagnate in an -neighbourhood of the fixed point [34]. One could loosely say that the super-linear convergence of is stronger than the truncation by .

  2. if the iteration by is linearly convergent with contraction factor , the truncated iteration will still converge linearly, but finally stagnate in an -neighbourhood of [48]. Again, one could loosely say that iteration by and truncation by balance each other, thus resulting in a larger neighbourhood of stagnation.

We shall assume that the truncation level has thus been chosen according to the desired re-construction accuracy and taking into account the possible influence due to the convergence behaviour of the iterator .

2.2 Preliminaries and basic algebraic operations

As already pointed out, we assume in general that the set is a vector space. Our example of the space of tensors of interest was already introduced in xmpl-motiv. This is clearly a vector space, so we are able to add two such tensors, and multiply each by some real number, in other words we may form linear combinations. The additive neutral element — the zero tensor — will be denoted as . It is important that the compressed storage format chosen allows to perform the vector space operations without going into the full representation. This will be shown for some of the familiar tensor formats in tensor-rep. Furthermore, as as vector spaces, we can carry the canonical Euclidean inner product on to , which for is denoted by

This makes into a Euclidean or Hilbert space. We assume that the computation of this inner product is feasible when both and are in a compressed representation; again this will be demonstrated in tensor-rep for the low-rank representations we shall consider as examples.

On this space an additional structure is needed, namely a multiplication which will make it into a unital associative and commutative algebra. In our example this will be the Hadamard multiplication, which apparently goes back to Schur, which is the point- or component-wise multiplication on . It is the usual point-wise multiplication of functions, in this case the functions . For , define the Hadamard product as


It is a bi-linear operation (linear in each entry), and it is commutative. Thus this product makes into an associative and commutative algebra. The symbol will be used both for the abstract algebra product on , as well as for the Hadamard product on . As is any such algebra, for any it holds that .

Especially, as a function of , the product is a linear map on . This is the familiar canonical representation of an associative algebra as linear maps on itself, i.e. in this case in a commutative sub-algebra of the algebra of all linear maps with concatenation as product.

It easy to see that this Hadamard algebra has a unit element for the product,


a tensor with all entries equal to unity — this makes into a unital algebra. Observe that by defining for all the all-ones vectors , the Hadamard unit has rank one according to eq:CPI: . Obviously one has , the identity in .

Having defined a unit or neutral element for multiplication, we say that has a multiplicative inverse iff there is an element, denoted by , such that


Obviously not all elements have an inverse, e.g. the zero element is never invertible. In any case, not every has a Hadamard inverse, for this it is necessary that all entries are non-zero, and then . Note that and as it should be, and it is easily seen that , the inverse of in .

Elements of the algebra which can be written as a square are called positive, they form a convex cone ; note that obviously the zero tensor and the multiplicative unit are positive, and that if an invertible is positive, so is its inverse . As usual, this defines an order relation , which implies that for any positive one has . Further one may observe that .

A certain compatibility of the inner product and the algebra product is needed, as we require that


i.e. that the action of the algebra product is self-adjoint, or, in other words, that the maps are self adjoint. This condition is satisfied for the Hadamard algebra.

The inner product with the unit element defines, as function of , a positive linear functional , as for one has

It is a kind of “trace” or un-normalised state functional on the algebra, as in the Hadamard algebra. In particular


Conversely, if the algebra comes equipped with such an un-normalised positive state functional, the eq:state-inn-pr may be taken as the definition of an inner product, which automatically satisfies eq:self-adj-mult.

We shall assume that we can compute the algebraic operations — at least approximately — in compressed representation (this will be shown for the Hadamard algebra for the low-rank formats in tensor-rep), and that we can compute the multiplicative inverse also in compressed representation — again at least approximately.

As in any unital algebra, if for the element fails to be invertible — this means that is not invertible in — we shall say that is in the spectrum of ; here — as we are in finite dimensional spaces — it is an eigenvalue, i.e. there is an eigenvector such that . From this discussion it is immediate that in the Hadamard algebra case, if is an eigenvalue of or , there must be a multi-index such that , and hence in the real Hadamard algebra one has that always . Thus the spectrum of resp.  is . Observe that according to eq:self-adj-mult, each representing map from above is self-adjoint, which means that also in general all the spectra are real — .

Specifically, the Euclidean Hadamard algebra which we have constructed is obviously isomorphic — as a Euclidean algebra — to with the canonical Euclidean inner product when equipped with the point- or element-wise Hadamard product . Let us denote this isomorphism by ; it implies some ordering of the terms of each .

More enlightening and less obvious may be the unital algebra isomorphism with the commutative sub-algebra of diagonal matrices in the full matrix algebra with the usual matrix multiplication. Let us denote this isomorphism by . If is the vector containing all the elements of the tensor , then is the corresponding diagonal matrix, i.e. . Hence, by choosing the canonical Euclidean basis in and its image by as a basis in , the canonical representation is itself represented by the matrix .

Let be another tensor, with associated and , then from the definition of the isomorphism :


which contains the important equality


from where one sees that


which helps in understanding the following algorithms. Note that this shows that the algebra with the Hadamard product is already jointly diagonalised and is essentially in its Gel’fand representation [60], which is an abstract way of saying what was already stated above, namely that the spectrum of an element are exactly the individual terms or data stored in the tensor. If in some abstract unital associative and commutative algebra the multiplication is not the point-wise multiplication — e.g. think of convolution — then it is known that modulo some technicalities [60] it is isomorphic via the Gel’fand “diagonalisation” morphism to a function algebra, and all the algorithms to follow would deal with the spectrum of instead of with the “values of ”. These two notions only coincide for a function algebra.

2.3 Post-processing Tasks

Here the different post-processing tasks are explained together with how they will be computed. This may involve a number of auxiliary functions. The computation of these — again through truncated iteration Algorithm 1 — and the way in which some of the computations can be enhanced or accelerated is shown in auxil.

Finding the maximum or minimum

of is the first task we consider. Observe, that the element of maximum modulus of the tensor is also equal to the -norm .

Finding the maximum means finding the index