P2663R5
Interleaved complex values support in std::simd

Published Proposal,

This version:
http://wg21.link/P2663R5.pdf
Authors:
(Intel)
(Intel)
Audience:
LEWG
Project:
ISO/IEC 14882 Programming Languages — C++, ISO/IEC JTC1/SC22/WG21

Abstract

The paper proposes to support interleaved complex values in std::simd

1. Motivation

ISO/IEC 19570:2018 (1) introduced data-parallel types to the C++ Extensions for Parallelism TS. [P1928R6] is proposing to make that a part of C++ IS. Intel supports the concept of a standard interface to SIMD instruction capabilities and we have made extra suggestions for other APIs and facilities for std::simd in document [P2638R0]. One of the extra features we suggested was the ability to work with interleaved complex values. It is now becoming common for processor instruction sets to include native support for these operations (e.g., Intel AVX512-FP16, ARM Neon) and we think that allowing std::simd to directly access these features is advantageous.

This document gives more details about the motivation for supporting complex-values, the different storage and implementation alternatives, and some suggestions for a suitable API.

2. Background

Complex-valued mathematics is widely used in scientific computing and many types of engineering, with applications including physical simulations, wireless signal processing, media processing, and much more besides. C++ already supports complex-valued operations using the std::complex template available in the <complex> header file and C supports complex values by using the _Complex keyword. C++ supports floating point elements, including std::float16 (C++23), float, double or long double, though in practice compilers permit integer complex values too. Such is the importance of interleaved complex values that some mainstream processor vendors now provide native hardware support for complex-value SIMD instructions (e.g., Intel AVX512-FP16, ARM Neon v8.3).

Any complex value is represented as a pair of values, one for the real component and one for the imaginary component. In C and C++ a complex value is a pair of two values of the appropriate data type, allocated to a single unit containing 2 elements which are stored contiguously. This is illustrated in the figure below. This storage layout is used in other languages too, allowing for easy data interchange between software written in different languages, or with comparable user-defined types. This format may also be used for hardware accelerators or interfaces too:

…

When many complex values are stored together, such as in a C-style array, a C++ std::array or a C++ std::vector, then the real and imaginary elements are interleaved in memory:

…

Many applications, libraries, programming languages and interfaces between diverse software and hardware components will store data in this interleaved format, and it is treated as the default layout for blocks of complex-valued data. To improve the compute efficiency of this data format both Intel and ARM have introduced native hardware support to allow efficient SIMD operations to be performed in this data format without having to resort to reformatting the data into a more SIMD-amenable form. Where the underlying target hardware does not have native support for interleaved complex values it is straight-forward to synthesize a sequence of operations which give the same effect.

Given the popularity of interleaved complex-values as a data exchange and compute format we propose that std::simd should be able to directly represent this simd format in a form exemplified as:

std::simd<std::complex<float>, 5>
std::simd<std::complex<double>, 9>

In all these cases the fundamental element type will be a suitable complex value, and the individual components of each complex element will be worked on as a single unit, according to the rules of complex arithmetic where appropriate. This includes:

One of the aims of std::simd is to make it possible to write a generic code, which can use either a scalar type or a data-parallel type interchangeably. To this end, std::simd provides overloads of many of the standard functions which operate in data-parallel mode. The same should apply to a simd of complex values; it should be possible to write an algorithm which uses either std::complex<T> or simd<std::complex<T>> with only a type replacement, not with major API changes. To this end std::simd should include some member functions, such as real or imag, to mirror those found in std::complex.

Note that an alternative way to create parallel complex values is to allow std::complex to use simd elements, such as std::complex<simd<float>>. This storage layout is called separated storage and discussion of this format is outside the scope of this paper.

3. Overview of updates required to support std::complex simd elements

There are three ways in which std::simd needs to change to accommodate std::complex elements:

Each of these points will now be discussed in more detail.

3.1. std::simd base modifications

In its current definition, std::simd allows the element type to be any cv-unqualified arithmetic type other than bool. The first modification is to extend the set of allowable element types to include any std::complex<> specialisation for a cv-unqualified floating-point type.

Any simd operation or function which relies only on an element being an atomic container which can be arbitrarily moved or duplicated as a single unit will continue to operate as expected. For example, all the following will work on simd<complex<T>> without modification:

simd_mask’s of complex values do not need to have any special behaviour since the mask is simply a collection of boolean values, and the underlying element type makes no difference.

Constructors, including generator functions, which work on atomic units (e.g., broadcasting constructor) will work without modification. Note that the existing broadcast constructor will already permit a real-valued object to be used in a broadcast, since an individual value may be converted to a complex value with a zero imaginary:

auto x = simd<std::complex<float>>(3.14f);
// [3.14 + 0i, 3.14 + 0i, ...]

Numeric operations (binary, unary, and compound-assignment) are mostly covered by [P1928R6] remark:

A simd object initialized with the results of applying the indicated operator to lhs and rhs as a binary element-wise operation

In the case of complex values this has the practical effect of ensuring that operators will perform their complex operations without changing the definition of std::simd (e.g., operator* will perform complex-valued multiplication). Furthermore, [P1928R6] notes that constraints will be applied as necessary on operators and functions:

Constraints: Application of the indicated operator to objects of type T is well-formed.

This will ensure that operators which are not permitted for complex numbers (e.g., operator%, operator<<, operator&) will be removed from the overload set, as will functions which rely on these operators (e.g., min, max, clamp).

3.2. Additional complex methods for std::simd

It is desirable for std::simd to be source- or template-replaceable with respect to its underlying element type. For example, given a simple code fragment which works with std::complex:

auto my_fn(auto cmplx_value)
{
    std::complex<float> rotator (0, 1);
    return (cmplx_value * rotator).real();
}

In this example the function will work with a std::complex<float> input. Ideally it should also work with an input of type simd<complex<float>>. For this to happen, firstly the operator* must work with a scalar value, and secondly it should be possible for the real() method to be invoked on a simd<complex<T>> value.

std::complex provides a number of member functions:

The assignment operator and the compound-assignment operators are already provided by std::simd so no further action is required. The remaining member functions that are required are considered in the following sub-sections.

3.2.1. Constructors

Conversion and copy methods already exist in std::simd and can be used for complex SIMD values. An additional constructor is needed to match the std::complex constructor which can build a complex value from separate real and imaginary components: constructed from real and imaginary simd values.

requires is_complex_number<T>
constexpr basic_simd(const basic_simd<T::value_type, UAbi>& r = {},
                     const basic_simd<T::value_type, UAbi>& i = {})

This constructor is constrained so that it can only be used for SIMD values with complex elements (some suitable concept will be added to enforce this). Like its std::complex counterpart it can be used to build a complex value from real components only, with insertion of zero imaginaries. This constructor can be used as in the following example:

simd<float, 8> reals = ...;
simd<float, 8> imaginaries = ...;

// Create a real-valued complex number (i.e., zero imaginaries)
simd<std::complex<float>, 8> realAsComplex (reals);

// Create a complex value from real and imaginary simd inputs
simd<std::complex<float>, 8> cmplx (reals, imaginaries);

3.2.2. Real/imag accessors

It should be possible to separately read or write the real and imaginary components. For example, given a simd of std::complex<float> values, the real components of the simd will be extracted as a value of type std::simd<float> instead. Similarly, given a simd<float> value those values can be inserted into the simd as the real components of each respective value.

The read accessors for std::complex come in two variants: a free function, and a member function. These would be defined in std::simd as follows:

// Free functions for accessing real/imaginary components
template <typename T, typename ABI> constexpr auto real(const basic_simd<complex<T>, ABI>&);
template <typename T, typename ABI> constexpr auto imag(const basic_simd<complex<T>, ABI>&);

// Member functions of std::simd for accessing real/imaginary components.
constexpr simd<T, size()> std::simd<std::complex<T>>::real() const;
constexpr simd<T, size()> std::simd<std::complex<T>>::imag() const;

The write accessors for std::complex are provided as member functions which allow the real and imaginary components of an existing complex value to be overwritten.

template<typename T, typename ABI>
constexpr void real(const basic_simd<T, ABI>& reals);

template<typename T, typename ABI>
constexpr void imag(const basic_simd<T, ABI>& reals);

The advantage of making real and imag members of std::simd is an API consistency between std::complex<T> and std::simd<std::complex<T>> for both getters and setters.

3.3. Additional free functions

The free functions for std::complex include numeric operators, such as +, -, * and /. These are already proposed in [P1928R6] and need not be considered further.

The remaining free functions which should be specialized for std::simd<complex<T>> include:

real(simd<complex>) returns the real part
imag(simd<complex>) returns the imaginary part
abs(simd<complex>) returns the magnitude of a complex number
arg(simd<complex>) returns the phase angle
norm(simd<complex>) returns the squared magnitude
conj(simd<complex>) returns the complex conjugate
proj(simd<complex>) returns the projection onto the Riemann sphere
polar(simd<complex>) constructs a complex number from magnitude and phase angle
exp(simd<complex>) complex base e exponential
log(simd<complex>) complex natural logarithm with the branch cuts along the negative real axis
log10(simd<complex>) complex common logarithm with the branch cuts along the negative real axis
pow(simd<complex>) complex power, one or both arguments may be a complex number
sqrt(simd<complex>) complex square root in the range of the right half-plane
sin(simd<complex>) computes sine of a complex number
cos(simd<complex>) computes cosine of a complex number
tan(simd<complex>) computes tangent of a complex numbe
asin(simd<complex>) computes arc sine of a complex number
acos(simd<complex>) computes arc cosine of a complex number
atan(simd<complex>) computes arc tangent of a complex number
sinh(simd<complex>) computes hyperbolic sine of a complex number
cosh(simd<complex>) computes hyperbolic cosine of a complex number
tanh(simd<complex>) computes hyperbolic tangent of a complex number
asinh(simd<complex>) computes area hyperbolic sine of a complex number
acosh(simd<complex>) computes area hyperbolic cosine of a complex number
atanh(simd<complex>) computes area hyperbolic tangent of a complex number

All of these functions should operate element-wise in the same way as their scalar counterparts. Note that although most of these functions take complex inputs and generate complex outputs, some of these functions generate real-valued outputs (e.g., abs returns the magnitude, which is real-valued), or accept an input which is real-valued (e.g., conj can accept a real-valued value and return a complex value).

4. Implementation experience

Intel have written an example implementation of std::simd which includes the extensions to support interleaved complex values. We have been able to generate efficient code both on targets with native support (e.g., AVX512_FP16), and also on targets without native support which require synthesized code sequences.

5. Wording

The wording relies on [P1928R6] being landed to the Working Draft. Below, substitute the � character with a number the editor finds appropriate for the table, paragraph, section or sub-section.

5.1. Modify [simd.general]

�: The set of vectorizable types comprises all cv-unqualified arithmetic types other than bool , and std::complex<> specializations of cv-unqualified floating-point types. .

5.2. Modify [simd.summary]

�: Default intialization performs no initialization of the elements , except when the value type is a std::complex<> specialization in which case all elements will be value-initialized ; value-initialization initializes each element with T().

5.3. Modify [simd.ctor]

simd constructors [simd.ctor]
[...]
template<typename CAbi>
constexpr explicit
basic_simd(const basic_simd<complex-element-type, CAbi> reals = {},
           const basic_simd<complex-element-type, CAbi> imags = {});

Constraints:

  • The value_type is a complex value (i.e., some specialization of std::complex<> for a cv-unqualified floating-point type).

  • The complex-element-type is a floating-point type which is the same as T::value_type (where T is the simd’s own value_type).

  • basic_simd<complex-element-type, CAbi>::size is equal to basic_simd::size.

Effects:

  • Construct a simd of elements of type value_type (where value_type is a specialization of std::complex<>), such that the ith element is equal to the object constructed as value_type(reals[i], imags[i]) for all i in the range [0..size).

Remarks:

  • Where no parameter values are provided this function acts as default initializer which value-initializes all real and imaginary elements.

5.4. Add new section [simd.complex_accessors]

simd complex accessors [simd.complex_accessors]
constexpr simd<complex-element-type, size> real() const;
constexpr simd<complex-element-type, size> imag() const;

Constraints:

  • T is a complex value (i.e., any valid instantiation of std::complex<> for a cv-unqualified floating-point type).

  • The complex-element-type is a floating-point type which is the same as T::value_type (where T is the simd’s own value_type).

Returns:

  • A simd with values of complex-element-type that is a value_type of std::complex of the same width. The ith element will be equal to (*this)[i].real() or (*this)[i].imag() respectively for all i in the range [0..size).

template<typename CAbi>
constexpr void real(const basic_simd<complex-element-type, CAbi>& v) noexcept;
template<typename CAbi>
constexpr void imag(const basic_simd<complex-element-type, CAbi>& v) noexcept;

Constraints:

  • T is a complex value (i.e., any valid instantiation of std::complex<> for a cv-unqualified floating-point type).

  • The complex-element-type is a floating-point type which is the same as T::value_type (where T is the simd’s own value_type).

  • basic_simd<complex-element-type, CAbi>::size is equal to basic_simd::size.

Effect:

  • Modify each element of the simd such that the real or imaginary component of the ith element is replaced by the value of v[i].real or v[i].imag() respectively.

template<typename T, typename Abi>
constexpr simd<complex-element-type, size>
real(const basic_simd<typename T, typename Abi>& v) noexcept;

template<typename T, typename Abi>
constexpr simd<complex-element-type, size>
imag(const basic_simd<typename T, typename Abi>& v) noexcept;

Constraints:

  • T is a complex value (i.e., any valid instantiation of std::complex<> for a cv-unqualified floating-point type).

  • The complex-element-type is a floating-point type which is the same as T::value_type (where T is the simd’s own value_type).

Returns:

  • A simd with values of complex-element-type that is a value_type of std::complex of the same width, where the ith element is equal to v[i].real() or v[i].imag() respectively for all i in the range [0..v.size).

5.5. Add new section [simd.complex_math]

simd complex maths library [simd.complex_library]

For each set of overloaded functions within <complex> there shall be additional overloads sufficient to ensure that:

  • Any basic_simd<std::complex<T>, Abi> type can be used where std::complex<> would be used (parameter or return value).

  • Any real-valued simd<T, Abi> type can be used where T would be used (parameter or return value).

It is unspecified whether a call to these overloads with arguments that are all convertible to basic_simd<T, Abi> but are not of type basic_simd<T, Abi> is well-formed.

Each function overload produced by the above rules applies the indicated <complex> function element-wise. For the mathematical functions, the results per element only need to be approximately equal to the application of the function which is overloaded for the element type.

The result is unspecified if a domain, pole, or range error occurs when the input argument(s) are applied to the indicated <complex> function. [ Note: Implementations are encouraged to follow the C specification (especially Annex F). — end note ]

6. Polls

6.1. Kona 2023 LEWG polls

Poll: Forward P2663R4 (Interleaved complex support in std::simd) with modified wording for 5.1 to constrain complex<> to floating point to LWG for C++26 (to be confirmed with a Library Evolution electronic poll).

SF F N A SA
8 8 0 0 0

6.2. Kona 2022 SG1 polls

Poll: After significant experience with the TS, we recommend that the next version (the TS version with improvements) of std::simd target the IS (C++26)

SF F N A SA
10 8 0 0 0

Poll: Future papers and future revisions of existing papers that target std::simd should go directly to LEWG. (We do not believe there are SG1 issues with std::simd today.)

SF F N A SA
9 8 0 0 0

7. Revision History

R4 => R5

R3 => R4

R2 => R3

R1 => R2

R0 => R1

References

Informative References

[P1928R6]
Matthias Kretz. std::simd - Merge data-parallel types from the Parallelism TS 2. 19 June 2023. URL: https://wg21.link/p1928r6
[P2638R0]
Daniel Towner. Intel's response to P1915R0 for std::simd parallelism in TS 2. 22 September 2022. URL: https://wg21.link/p2638r0