- Generated by
1.12.0
|
TRIQS/nda 1.3.0
Multi-dimensional array library for C++
|
In this example, we give a quick overview of the basic capabilities of nda.
All the following code snippets are part of the same main function:
In its simplest form, an nda::array has 2 template parameters
The following creates an integer (int) array of rank 2 and initializes its elements:
The first line calls the constructor with the intended shape of the array, i.e. the extents along the two dimensions. In this case, we want it to be 3-by-2.
In the second line, we use a range-based for loop (note the reference auto &) to assign the first few positive integer numbers to the elements of the array.
We could have achieved the same using the constructor which takes std::initializer_list objects:
By default, nda::array stores its elements in C-order. To create an array in Fortran-order, we can specify a third template parameter:
Here, nda::F_layout is one of the Layout policies.
While in 2-dimensions, the only possibilities are C-order or Fortran-order, in higher dimensions one can also specify other stride orders (see nda::basic_layout and nda::basic_layout_str).
Let's check the contents, sizes and shapes of the arrays using the overloaded streaming operators:
Output:
You can see the difference between the memory layouts of the array:
Note: The reason why we see a difference in the two arrays is because the range-based for loop is optimized to iterate over contiguous memory whenever possible. If we would have used traditional for loops instead,
int n = 0;for (int i = 0; i < 3; ++i) {for (int j = 0; j < 2; ++j) {B(i, j) = n++;}}there would be no difference between the output of
AandB:B =[[0,1][2,3][4,5]]
We can access single elements of the array using the function call operator of the array object. For a 2-dimensional array, we have to pass exactly two indices, otherwise it won't compile:
Output:
Note: Out-of-bounds checks are only enabled in debug mode. The following code will compile and result in undefined behavior:
std::cout << "A(3, 2) = " << A(3, 2) << std::endl;Output (might be different on other machines):
A(3, 2) = 9
It is straightforward to assign a scalar or another array to an existing array:
Output:
Views offer a lightweight and efficient way to manipulate and operate on existing arrays since they do not own their data, i.e. there is no memory allocation or copying involved when creating a view (see nda::basic_array_view).
Taking a view on an existing array can be done with an empty function call:
Output:
Here, A_v has the same data, size and shape as A. Since A_v is a view, it points to the same memory location as A. That means if we change something in A_v, we will also see the changes in A:
Output:
In most cases, views will just behave like arrays and the majority of functions and operations that can be performed with arrays, also work with views.
A slice is a view on only some parts of an existing array.
To take a slice of an array, one uses again the function call operator but this time with one or more nda::range, nda::range::all or nda::ellipsis objects:
nda::range is similar to a Python range and specifies which elements should be included in the slice along a certain dimension.nda::range::all is similar to Python's : syntax and specifies that all elements along a certain dimension should be included in the slice.... syntax and specifies that all "missing" dimensions should be included in the slice. One can achieve the same with multiple nda::range::all objects.For example, we can view only the second row of A:
Output:
As seen in this example, slices usually have a different size and shape compared to the original array. In this case, it is a 1-dimensional view of size 2.
Since a slice is just a view, everything mentioned above about views still applies:
Output:
We can perform various Arithmetic operations with arrays and views. Arithmetic operations are (mostly) implemented as lazy expressions. That means the operations are not performed right away but only when needed at a later time.
For example, the following tries to add two 2-dimensional arrays element-wise:
Output:
Note that C + D does not return the result of the addition in the form of a new array object, but instead it returns an nda::expr object that describes the operation it is about to do.
To turn a lazy expression into an array/view, we can either assign it to an existing array/view, construct a new array or use nda::make_regular:
Output:
The behavior of arithmetic operations (and some other functions) depends on the algebra (see nda::get_algebra) of its operands.
For example, consider matrix-matrix multiplication ('M' algebra) vs. element-wise multiplication of two arrays ('A' algebra):
Output:
Here, an nda::matrix is the same as an nda::array of rank 2, except that it belongs to the 'M' algebra instead of the 'A' algebra.
Similar to arithmetic operations, most of the Mathematical functions that can operate on arrays and views return a lazy function call expression (nda::expr_call) and are sensitive to their algebras:
Output:
Besides arithmetic operations and mathematical functions, one can apply some basic Algorithms to arrays and views:
Output:
Most of the algorithms in the above example are self-explanatory, except for the nda::any and nda::all calls. They expect the given array elements to have a bool value type which is not the case for the C array. We therefore use nda::map to create a lazy nda::expr_call that returns a bool upon calling it. Since nda::expr_call fulfills the nda::Array concept, it can be passed to the algorithms.
Writing arrays and views to HDF5 files using nda's HDF5 support is as easy as
Dumping the resulting HDF5 file gives
Reading the same file into a new array is just as simple:
Output:
nda has a BLAS interface and an LAPACK interface and provides the nda::matrix and nda::vector types. While nda::matrix is a 2-dimensional array belonging to the 'M' algebra, nda::vector is a 1-dimensional array belonging to the 'V' algebra. Under the hood, both of them are just type aliases of the general nda::basic_array class. Everything mentioned above for nda::array objects is therefore also true for matrices and vectors, e.g. writing/reading to HDF5, accessing single elements, taking slices, etc.
Let us demonstrate some of the supported Linear algebra features:
Output:
Note: Matrix-matrix and matrix-vector multiplication do not return lazy expressions, since they call the corresponding BLAS routines directly, while element-wise array-array multiplication does return a lazy expression.
nda contains the Compile-time lazy expressions and functions (CLEF) library which is a more or less standalone implementation of general lazy expressions.
One of the most useful tools of CLEF is its Automatic assignment feature:
Output:
Here, we first construct an uninitialized 4-by-7 array before we use automatic assignment to initialize the array. This is done with the overloaded operator<< and nda::clef::placeholder objects (i_ and j_ in the above example). It will simply loop over all possible multidimensional indices of the array, substitute them for the placeholders and assign the result to the corresponding array element, e.g. F(3, 4) = 3 * 7 + 4 = 25.
This is especially helpful for high-dimensional arrays where the element at (i, j, ..., k) can be written as some function \( g \) of its indices, i.e. \( F_{ij \dots k} = g(i, j, \dots, k) \).
The above features constitute only a fraction of what you can do with nda.
For more details, please look at the other Examples and the API Documentation.
1.12.0