In-depth: Tensors, ndarrays for .NET

Version 9 of Numerics.NET has just been released. The biggest new feature in this version is the addition of a tensor library.

Tensors or multi-dimensional arrays are a fundamental building block of numerical computing.

Python has long had the numpy library. For years it has been at the core of Python’s scientific computing ecosystem. Libraries like pandas and scikit-learn are largely built on top of numpy. The main object in numpy, and the main source of its power, is a multi-dimensional array.

With this release, Numerics.NET brings much of the functionality of numpy to the .NET platform. This includes:

  • A rich multi-dimensional array (tensor) type.
  • Powerful features like broadcasting and advanced indexing.
  • A comprehensive set of mathematical functions, linear algebra, and more.
  • Easy interop with both managed arrays and external or native arrays.
  • Excellent performance with the infrastructure in place to grow even faster.

While numpy has obviously been a source of inspiration, we deliberately chose not to simply copy the API but do what makes most sense on the .NET platform.

For example, the tensor type is generic over its element type. This follows the pattern of just about every other .NET collection type since the introduction of generics in .NET 2.0.

We’ve also followed a pattern we’ve used elsewhere in Numerics.NET: The main type of the tensor library is the Tensor<T> class. There is also a static class, Tensor, that contains methods to create tensors and perform operations on them.

Getting started

The tensor library has its own NuGet package, Numerics.NET.Tensors, but is included as a reference in the main Numerics.NET package. The main classes live in the Numerics.NET namespace, while classes that support the tensor library live in the Numerics.NET.Tensors namespace.

The tensor library is designed to be easy to use. Here is a simple example:

C# 
using Numerics.NET;

// Create a tensor (1 x 4).
var t0 = Tensor.Create([1, 2, 3, 4]); // [1, 2, 3, 4]

// Reshape tensor (4 x 1).
var t1 = t0.Reshape(4, 1);            // [[1], [2], [3], [4]]

Indexers always return tensors. To get or set an individual element, use the GetValue and SetValue methods:

C# 
// Slice tensor (2 x 1).
var t2 = t1[1.., ..]; // [[2], [3]]
// Slice tensor (2 x 1).
var t3 = t0[1]; // [2] (Tensor<int>)

// Get a single value.
var s1 = t0.GetValue(1); // 2 (int)
var s2 = t1.GetValue(0, 1); // 2 (int)
var s3 = t1.GetValue(0); // Runtime error: RankException

Indexers return a view of the original tensor, not a copy. The exception is when sets of integers or a boolean mask is used as an index. This doesn’t apply when assigning to a tensor. For example, you can conditionally set values in a tensor:

C# 
t1[t1 > 2] = 99; // [[1], [2], [99], [99]]

Tensors support a wide range of operations. The standard arithmetic operators are overloaded:

C# 
// Math operations.
var t4 = t0 + 1; // [2, 3, 4]
var t5 = t0 + t0; // [2, 4, 6]
var t6 = t0 - 1; // [0, 1, 2]
var t7 = t0 - t0; // [0, 0, 0]
var t8 = t0 * 2; // [2, 4, 6]
var t9 = t0 * t0; // [1, 4, 9]
var t10 = t0 / 2; // [0.5, 1, 1.5]
var t11 = t0 / t0; // [1, 1, 1]

Broadcasting - expanding the singleton dimension of a tensor to match the shape of another tensor - is automatic.

C# 
// Broadcast (3 x 1) + (1 x 3) tensors -> (3 x 3).
var t12 = t0 + t1;
// [
//  [ 2, 3, 4],
//  [ 3, 4, 5],
//  [ 4, 5, 6]
// ]

The Tensor class has equivalent methods for each operator that offer additional options. For example, you can specify the destination of an operation, like Tensor.Add(left, right, result). This avoids unnecessary allocations and can dramatically improve performance. You can also provide a condition that must be satisfied for an operation to be applied.

C# 
// Specify result. Note that t5 is alread the correct shape.
Tensor.Add(t0, 1, t5); // t5 -> [2, 3, 4]
// Specify condition
var tb = t0 % 2 == 1; // [true, false, true]
Tensor.Add(t0, 10, t5, mask: tb); // t5 -> [11, 2, 13]

There are various ways of manipulating the shape of a tensor. All of these operations return a view of the original tensor, not a copy, whenever possible.

C# 
var t13 = Tensor.CreateRange(12).Reshape(3, 4);
// [[0, 1, 2, 3], [4, 5, 6, 7], [8, 9, 10, 11]]
var t14 = t13.Transpose(); 
// [[0, 4, 8], [1, 5, 9], [2, 6, 10], [3, 7, 11]]
var t15 = t14.InsertDimension(1);
// t15.Shape -> [4, 1, 3]
var t16 = t15.SwapAxes(1, 2);
// t16.Shape -> [4, 3, 1]
var t17 = t16.Squeeze(^1);
// t17.Shape -> [4, 3]

Some unique features

The tensor library has a number of features that are unique to Numerics.NET:

Grouped reductions

Reductions are great for things like summing or finding the maximum value of the elements of a tensor along a dimension. But what if you want to do something more complex, like computing the mean of a set of points that are grouped together?

A practical application is found in clustering, where we want to compute the centroid of each cluster. The coordinates of the centroid of a cluster are the means of the corresponding coordinates of the cluster members.

Let’s say we have a 2D tensor points where each row corresponds to a point. We also have an integer tensor clusters that specifies which cluster each point belongs to.

Traditionally, there are two ways to achieve this. The first is to use advanced indexing. The idea is to loop through the clusters, select the points that belong to each cluster, and then compute their mean to obtain the centroid.

C# 
var centroids = Tensor.Create<double>((nClusters, points.Shape[1]));
for (int i = 0; i < nClusters; i++)
{
    var mask = clusters == i;
    var clusterPoints = points[mask];
    centroids[i] = Tensor.Mean(clusterPoints, axis: 0);
}

Alternatively, we can use a conditional operation. We still loop through the clusters, but instead of creating a whole new tensor, we apply the operation to the original tensor with a mask that selects the points that belong to the current cluster.

C# 
var centroids = Tensor.Create<double>((nClusters, points.Shape[1]));
for (int i = 0; i < nClusters; i++)
{
    var mask = clusters == i;
    centroids[i] = Tensor.Mean(points, axis: 0, mask: mask);
}

Although the second alternative is an improvement over the first - there is much less memory allocation - it is still not ideal. We are still looping through the clusters and so we’re going through the whole set of points multiple times.

In Numerics.NET, the solution is much simpler. Reduction methods have an overload that takes a grouping parameter. This parameter specifies which elements belong to the same group. Instead of returning a single value along each axis, it returns the set of reduced values for each group:

C# 
var centroids = Tensor.Mean(points, axis: 0, groupBy: clusters);

Views

Tensors often contain millions of elements. Creating a new tensor with its own storage just to create a view of a subset of the original tensor is wasteful. Numerics.NET tries to avoid this as much as possible.

Nearly all structural operations on a tensor return a view of the original tensor and don’t copy the underlying data. When you slice a tensor, you get a view of the original tensor. When you transpose a tensor or move its axes around, you get a view of the original tensor. When you add or remove singleton dimensions, you get a view of the original tensor.

There are only a handful of exceptions to this rule:

  1. Advanced indexing (e.g. using a mask or a set of integers as an index), always copies the elements.
  2. When you reshape a tensor, it is not always possible to create a view.

In addition, Numerics.NET has a couple of mechanisms that further reduce the need for copying.

Lazy transformations

Lazy transformations are a way to apply an operation to an operand without explicitly forming the resulting tensor.

Mixed operations

Tensors in Numerics.NET are strongly typed. This means that you can’t add just add a tensor of integers to a tensor of doubles. However, you can cast the tensor of integers to a tensor of doubles and then add them.

In Numerics.NET, casting is done using the As<T>() method. This method returns a lightweight view of the original tensor with a different element type. The values aren’t copied. Instead, the operation is performed lazily.

C# 
var t18 = t0.As<double>() + 1.0;

The only allocation here, aside from the result tensor, is the memory needed to store the transformed tensor object.

Conditional operations

When you apply an operation to a tensor, you can specify a condition that must be satisfied. Some operations, particularly logical operations or relational operations involving scalars, are again performed lazily. When used as a mask in an operation or an advanced indexing operation, the condition is evaluated only when needed.

C# 
var t19 = Tensor.Multiply(t0, 2, mask: t0 > 2);

Complex tensors

Obtaining the real and imaginary parts of a complex tensor is another example of a lazy transformation that illustrates a further point: some transformations can go both ways.

The RealPart and ImaginaryPart methods return views of a complex tensor that allow you to work with the real and imaginary parts independently.

C# 
var tc = Tensor.CreateFromFunction<Complex<int>>(3, (i) => (i, 2 * i + 1));
// [ 0 + 1i, 1 + 3i, 2 + 5i ]
var tcReal = tc.RealPart();
// [ 0, 1, 2 ]
var tcImag = tc.ImaginaryPart();
// [ 1, 3, 5 ]

Not only is the tensor view updated when the original tensor is changed, but when you update a view, the original tensor is updated as well:

C# 
tc.SetValue((8, 9), 1);
Console.WriteLine($"{tc[1]} == ({tcReal[1]}, {tcImag[1]}");
// Likewise, when you update a view, the original tensor is updated as well:
tcReal.SetValue(6, 1);
// tc[1] = (6, 9)

You can even do fun things like swap the real and imaginary parts:

C# 
Tensor.Swap(tcReal, tcImag);
// tc -> [ 1 + 0i, 9 + 6i, 5 + 2i ]

Conclusion

The tensor library in Numerics.NET is a powerful tool for numerical computing. It is designed to be easy to use, efficient, and flexible.

It is also designed to be a foundation for future development. We have plans to add more features, including:

  • Support for sparse tensors.
  • Reductions over multiple axes.
  • Lazy transformations involving more than one argument.
  • Axis indexes and labels, similar to data frames.
  • Automatic differentiation.
  • Many more performance optimizations.

The infrastructure for these features is already in place. We just didn’t get to implementing them yet.

If you have ideas or suggestions, please let us know. We are always looking for ways to improve Numerics.NET and make it more useful to you.