An Overview Of C++26: The Library – Math

An Overview Of C++26: The Library – Math

std::submdspan

std::submdspan is s a subspan of an existing span std::mdspan (C++23). The subspan didn’t make it into C++23 and was added with C++26.

Before I continue with C++26, I must make a short detour to C++23. hi

std::mdspan

A std::mdspan is a non-owning multidimensional view of a contiguous sequence of objects. The contiguous sequence of objects can be a plain C-array, a pointer with a size, a std::array, a std::vector, or a std::string. Often, this multidimensional view is called a multidimensional array.

The number of dimensions and the size of each dimension determine the shape of the multidimensional array. The number of dimensions is called rank, and the size of each dimension extension. The size of the std::mdspan is the product of all dimensions that are not 0. You can access the elements of a std::mdspan using the multidimensional index operator [].

Each dimension of a std::mdspan can have a static or dynamic extent. static extent means that its length is specified at compile time; dynamic extent means that its length is specified at run time.

Thanks to class template argument deduction (CTAG) in C++17, the compiler can often automatically deduce the template arguments from the types of initializers:

// mdspan.cpp

#include <mdspan>
#include <iostream>
#include <vector>

int main() {
    
    std::vector myVec{1, 2, 3, 4, 5, 6, 7, 8};          // (1)

    std::mdspan m{myVec.data(), 2, 4};                  // (2)
    std::cout << "m.rank(): " << m.rank() << '\n';      // (4)

    for (std::size_t i = 0; i < m.extent(0); ++i) {     // (6)
        for (std::size_t j = 0; j < m.extent(1); ++j) { // (7)
            std::cout << m[i, j] << ' ';                // (8)
        }
        std::cout << '\n';
    }

    std::cout << '\n';

    std::mdspan m2{myVec.data(), 4, 2};                 // (3)
    std::cout << "m2.rank(): " << m2.rank() << '\n';    // (5)

    for (std::size_t i = 0; i < m2.extent(0); ++i) {
        for (std::size_t j = 0; j < m2.extent(1); ++j) {
        std::cout << m2[i, j] << ' ';  
    }
    std::cout << '\n';
  }

}
        

I apply class template argument deduction three times in this example. Line (1) uses it for a std::vector, and lines (2) and (3) for a std::mdspan. The first 2-dimensional array m has a shape of (2, 4), the second one m2 a shape of (4, 2). Lines (4) and (5) display the ranks of both std::mdspan. Thanks to the extent of each dimension (lines 6 and 7) and the index operator in line (8), it is straightforward to iterate through multidimensional arrays.


Here ends my detour, and I continue with std::submdspan.

std::submdspan

This function was considered critical for the overall functionality of mdspan. However, due to review time constraints, it was removed in order for mdspan to be included in C++23.

?


Modernes C++ Mentoring

Do you want to stay informed: Subscribe.

?

Creating a std::submdspan is straightforward. Its first parameter is an mdspan x, and the remaining x.rank() parameters are slice specifiers, one for each dimension of x. The slice specifiers describe which elements of the range [0,x.extent(d)) are part of the multidimensional index space of the returned mdspan.

This leads to the following fundamental signature:

template<class T, class E, class L, class A,
         class ... SliceArgs)
auto submdspan(mdspan<T,E,L,A> x, SliceArgs ... args);        

where E.rank() must be equal to sizeof...(SliceArgs).

Proposal P2630R4 provides, in addition to the definition of a std::submdspan, a few examples for a rank-1 mdspan.

int* ptr = ...;
int N = ...;
mdspan a(ptr, N);

// subspan of a single element
auto a_sub1 = submdspan(a, 1);
static_assert(decltype(a_sub1)::rank() == 0);
assert(&a_sub1() == &a(1));

// subrange
auto a_sub2 = submdspan(a, tuple{1, 4});
static_assert(decltype(a_sub2)::rank() == 1);
assert(&a_sub2(0) == &a(1));
assert(a_sub2.extent(0) == 3);

// subrange with stride
auto a_sub3 = submdspan(a, strided_slice{1, 7, 2});
static_assert(decltype(a_sub3)::rank() == 1);
assert(&a_sub3(0) == &a(1));
assert(&a_sub3(3) == &a(7));
assert(a_sub3.extent(0) == 4);

// full range
auto a_sub4 = submdspan(a, full_extent);
static_assert(decltype(a_sub4)::rank() == 1);
assert(a_sub4(0) == a(0));
assert(a_sub4.extent(0) == a.extent(0));
        

The same rules apply to the multi-dimensional use case.

int* ptr = ...;
int N0 = ..., N1 = ..., N2 = ..., N3 = ..., N4 = ...;
mdspan a(ptr, N0, N1, N2, N3, N4);

auto a_sub = submdspan(a,full_extent_t(), 3, strided_slice{2,N2-5, 2}, 4, tuple{3, N5-5});

// two integral specifiers so the rank is reduced by 2
static_assert(decltype(a_sub) == 3);
// 1st dimension is taking the whole extent
assert(a_sub.extent(0) == a.extent(0));
// the new 2nd dimension corresponds to the old 3rd dimension
assert(a_sub.extent(1) == (a.extent(2) - 5)/2);
assert(a_sub.stride(1) == a.stride(2)*2);
// the new 3rd dimension corresponds to the old 5th dimension
assert(a_sub.extent(2) == a.extent(4)-8);

assert(&a_sub(1,5,7) == &a(1, 3, 2+5*2, 4, 3+7));
        

This is not the end of the support is C++26.

<linalg>:

<linalg> is a free function linear algebra interface based on the BLAS.

BLAS: Basic Linear Algebra Subprograms (BLAS) is a specification that prescribes a set of low-level routines for performing common linear algebra operations such as vector addition, scalar multiplication, dot products, linear combinations, and matrix multiplication. They are the de facto standard low-level routines for linear algebra libraries … .(https://en.wikipedia.org/wiki/Basic_Linear_Algebra_Subprograms)

Proposal P1673R13 proposes a C++ Standard Library dense linear algebra interface based on the dense Basic Linear Algebra Subroutines (BLAS). This corresponds to a subset of the BLAS Standard.

This linear algebra was long missed in the C++ community. According to the proposal the following code snippet is the hello world of linear algebra. It scales the elements of a 1-D mdspan by a constant factor, first sequentially, then in parallel.

  constexpr size_t N = 40;
  std::vector<double> x_vec(N);

  mdspan x(x_vec.data(), N);
  for(size_t i = 0; i < N; ++i) {
    x[i] = double(i);
  }

  linalg::scale(2.0, x); // x = 2.0 * x
  linalg::scale(std::execution::par_unseq, 3.0, x);
  for(size_t i = 0; i < N; ++i) {
    assert(x[i] == 6.0 * double(i));
  }
        

What’s next?

I’m still not done with C++26. In my next post, I will write about saturation arithmetic and the concurrency support of C++26.

要查看或添加评论,请登录

社区洞察

其他会员也浏览了