The Tensor class

Basics of the Tensor class

A Tensor is a multidimensional array of numbers; in practice, they are constrained to real or complex double precision numbers. Their behavior is similar to Matlab's arrays in that they can store only numbers, be accessed with one or more indices using the () or [] syntaxes, reshaped, sliced, and all that with an automated memory management.

Internally, a Tensor consists of two data items: A Vector containing the dimensions, and a Vector containing the raw data. The Tensor class (or rather the underlying Vector class) makes excessive use of copy-on-write. Therefore, modifying the dimensions of a tensor is usually a cheap operation, while modifying the data may be expensive.

To simplify the interfacing to BLAS/LAPACK libraries, the data is aligned in row-major ordering (as in Fortran, and in contrast to C++).

For all functions that require indices of the tensor (e.g., accessing elements, creating a tensor), explicit functions for up to six-dimensional tensors are provided. For cases with more than six dimensions, or if the dimensions should not be required, there is an additional function that works with an Indices class. The setup is explained further below.

Creating a Tensor

To create a new tensor with uninitialized content, you can use the standard constructor

CTensor ca = CTensor(2, 3);    // creates a 2x3 matrix of complex doubles
RTensor ra = RTensor(5, 4, 3); // creates a 5x4x3 tensor of doubles

For convenience, there exist a couple of other static functions that produce standard tensors

RTensor a = RTensor::random(3, 2);  // fill the tensor with random content
CTensor b = CTensor::eye(5, 5);     // unity matrix; only works for up to 2 dimensions
CTensor c = CTensor::zeros(3, 4, 5);// fills the tensor with zeros
RTensor d = RTensor::ones(4, 3);    // fills the tensor with ones; explicit functions only for up to 2 dims

Tensors can also be created statically inside the code; this is explained further below. Finally, there is also some basic functionality to write Tensors to files and read them afterwards; for this, see the sdf namespace.

Accessing single elements of a tensor

There are two different access mechanisms if you want to retrieve a single entry in the Tensor. You can view the data as a linear sequence of numbers and read a given index, or you can use the dimension information of the tensor. For each access mechanism, a reading and a writing function is provided.

To access the data sequentially, you can use the square brackets, and the function at_seq().

RTensor t = RTensor::eye(2,2);  // The final tensor is ( 1 0 )
t.at(1,0) = 2;                  //                     ( 2 1 )
t[0];               // returns 1.0
t[1];               // returns 2.0
t.at_seq(2) = 5.1;  // sets the upper right number

To access the data using the Tensor's dimensions, you can use brackets or the function at(). Note that indices can be positive or negative, where the latter means "count from the end". For example, -1 means the last element along a dimension. This convention is by the way also used whenever a dimension is requested, a -1 would refer to the last dimension.

RTensor t = RTensor::eye(2,2);
t(0,0);                        // first row, first column
t(1,0);                        // second row, first column
t(-1,1);                       // last (i.e., second) row, second column
t.at(1,0) = t.at(0,1) = 1;     // now all elements of t are 1.

The various write functions return a reference to the corresponding element. To disallow these functions, you can always declare a tensor as constant.

const RTensor t = RTensor::random(5, 4);
t[0];                                    // OK
t(4,0);                                  // OK
t.at_seq(3) = 1.2;                       // compilation error: modifies tensor
t.at(2,2)   = 5.0;                       // also fails on compilation

Tensor slicing

To access more than one element of a tensor, you can use the range() function together with ordinary brackets to get slices. This works similar to Matlab's slice notation. The range() function can take a number of integer arguments, or an Indices vector (see the section on using Indices) to access various items. As with single entries, negative values are accepted, and count from the end of the dimension. Note that all parameters either have to be numbers, thus retrieving a single item, or slices.

CTensor t = CTensor::random(5,5);

CTensor full  = t(range(), range());          // no argument  -> take all Indices
CTensor row   = t(range(2), range());         // one argument -> take only given index
CTensor full2 = t(range(1,-1), range());      // two arguments-> start and stop
CTensor evenRows = t(range(1,-1,2), range()); // three args   -> start,stop,stride
CTensor oddRows = t(range(igen<<1<<3), range());// Indices as argument

To assign data, you can use again the at() function. You can assign three things: Either a single value or another tensor or slice with the same dimensions.

CTensor dest = CTensor::zeros(5,3);
CTensor col  = CTensor::ones(5,1);

dest.at(range(), range(0)) = 5.0;     // first column gets complex value (5.0, 0)
dest.at(range(), range(1)) = col;     // second column gets content of col
dest.at(range(), range(2)) = dest(range(), range(1));  // copy to third column

Note that the return value is in both cases an internal data structure that is either transparently cast to a tensor, or accepts the assigned data. If you keep this data without casting to another Tensor, you will effectively create a loophole in the copy-on-write mechanism. For this reason, you must never use the C++-11 auto feature, unless you know exactly what you are doing (you do not).

CTensor t = CTensor::random(5,5);
auto view = t(range(), range(0,1));    // creates an unsafe object that bypasses copy-on-write

Tensor shapes and ranks

To query the dimensions, the tensor class offers a couple of functions:

CTensor t = CTensor::zeros(5, 4, 3);

index rank = t.rank();           // number of dimensions: 3
index size = t.size();           // total size: 5*4*3 = 60

int d1 = t::dimension(0);        // number of entries in first dimension: 5
int d2 = t::dimension(1);        // along second dimension: 4
int d3 = t::dimension(2);        // along third dimension: 3

t::get_dimensions(&d1, &d2, &d3);// same: returns the three dimensions

Indices dims = t::dimensions();  // Returns an index vector with the dimensions
d1 = dims[0];                    // etc.

The function dimensions() is handy to produce an Indices vector that can later be used for convenient manipulations if you want to abstract away the exact rank of the tensor. For two-dimensional tensors, and only for those, the class also offers the rows() and columns() functions, which are equivalent to dimension(0) and dimension(1);

If you want to modify the shape of a tensor, you can use the reshape() function

CTensor t = CTensor::random(5,4);
CTensor t2 = reshape(t, 4, 5);

t.dimension(0) == t2.dimension(1);  // both comparisons evaluate to true
t.dimension(1) == t2.dimension(0);

Using Indices vectors

Typically, the library provides functions that are overloaded for up to six dimensional tensors. Sometimes, however, you might have larger objects, or you want to write a function or routine that can deal with tensors of varying dimensions. In this case, you have to work with Indices.

At its heart, an Indices object is just a fixed-sized vector of integer values. It uses the same copy-on-write mechanism as the Tensor class (actually, the data of tensors is implemented with the same vector class). You can use again the square brackets and the at() function with an index parameter for read / write access to the content.

CTensor t = CTensor::random(4, 5, 6);
Indices dims = t.dimensions();         // dims and t.dimensions() share the data

dims.at(0) = dims[2];                  // copy on write: now dims has its own copy
dims.at(2) = 4;

CTensor t2 = reshape(t, dims);

There are various ways to construct an Indices object. Either you get it from somewhere else (like Tensor::dimensions()), or you create it by claiming some memory or static initialization.

Indices i1 = Indices(5);              // allocate 5 entries for further use
Indices i2 = igen << 2 << 2 << 4      // creates an Indices object with entries (2,2,4); see next section
Indices i3 = Indices::range(1, 5, 2)  // start,stop,stride as parameters; i3 = (1, 3, 5)
Indices i4 = i2 << i3;                // concatenation; i4 = (2, 2, 4, 1, 3, 5)

The functions all_equal() and some_unequal() can be used to compare Indices. Furthermore, all comparison operators (==, <= etc.) are overloaded, and will return a vector of booleans that gives the elementwise result of the comparison.

Indices i1 = igen << 1 << 2 << 3;
Indices i2 = igen << 1 << 2 << 5;

bool eq = all_equal(i1, i2);      // evaluates to false
Boolean comp = (i1 == i2)         // comp = (true, true, false)

Statically creating tensors

It is possible to create a tensor or vector using compile-time expressions. The syntax always follows the same pattern: first the generator, then the content separated by the operator <<:

RTensor t1 = rgen << 1.0 << 2.0 << 3.0;   // real-valued vector with elements (1,2,3)
CTensor t2 = cgen << cdouble(1.0, 2.0);   // complex-valued vector with element (1+2i)
Indices i1 = igen << 1 << 2;              // Indices object with entries (1,2)
Booleans b = bgen << true << false;       // boolean vector with elements (true,false)

A special generator is "xgen", which will convert automatically to the type of the first object that is fed into it. Internally, these expressions are evaluated using recursive templates; consequently they are evaluated at compile time and will slow down compilation if used excessively (e.g., if you feed 1000 elements in this way).

Note that the direct use of the generators always results in one-dimensional tensors of appropriate size. If you want to generate a tensor of a given size, you have to supply the dimensions as Indices vector. Typically, you will create this vector also statically. Note that the data is as usual interpreted in row-major form.

RTensor a(rgen << 1 << 2 << 3 << 4,   // creates matrix ( 1   2 )
          igen << 2 << 2);            //                ( 3   4 )