Reference documentation for deal.II version 9.5.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
Classes | Namespaces | Functions
vectorization.h File Reference
#include <deal.II/base/config.h>
#include <deal.II/base/exceptions.h>
#include <deal.II/base/template_constraints.h>
#include <array>
#include <cmath>

Go to the source code of this file.

Classes

struct  EnableIfScalar< VectorizedArray< Number, width > >
 
class  VectorizedArrayIterator< T >
 
class  VectorizedArrayBase< T, width >
 
class  VectorizedArray< Number, width >
 
struct  internal::VectorizedArrayTrait< T >
 
struct  internal::VectorizedArrayTrait< VectorizedArray< T, width_ > >
 
struct  std::iterator_traits<::VectorizedArrayIterator< T > >
 

Namespaces

namespace  internal
 
namespace  std
 STL namespace.
 

Functions

template<typename Number , std::size_t width>
inline ::VectorizedArray< Number, width > std::sin (const ::VectorizedArray< Number, width > &x)
 
template<typename Number , std::size_t width>
inline ::VectorizedArray< Number, width > std::cos (const ::VectorizedArray< Number, width > &x)
 
template<typename Number , std::size_t width>
inline ::VectorizedArray< Number, width > std::tan (const ::VectorizedArray< Number, width > &x)
 
template<typename Number , std::size_t width>
inline ::VectorizedArray< Number, width > std::exp (const ::VectorizedArray< Number, width > &x)
 
template<typename Number , std::size_t width>
inline ::VectorizedArray< Number, width > std::log (const ::VectorizedArray< Number, width > &x)
 
template<typename Number , std::size_t width>
inline ::VectorizedArray< Number, width > std::sqrt (const ::VectorizedArray< Number, width > &x)
 
template<typename Number , std::size_t width>
inline ::VectorizedArray< Number, width > std::pow (const ::VectorizedArray< Number, width > &x, const Number p)
 
template<typename Number , std::size_t width>
inline ::VectorizedArray< Number, width > std::pow (const ::VectorizedArray< Number, width > &x, const ::VectorizedArray< Number, width > &p)
 
template<typename Number , std::size_t width>
inline ::VectorizedArray< Number, width > std::abs (const ::VectorizedArray< Number, width > &x)
 
template<typename Number , std::size_t width>
inline ::VectorizedArray< Number, width > std::max (const ::VectorizedArray< Number, width > &x, const ::VectorizedArray< Number, width > &y)
 
template<typename Number , std::size_t width>
inline ::VectorizedArray< Number, width > std::min (const ::VectorizedArray< Number, width > &x, const ::VectorizedArray< Number, width > &y)
 
Packing and unpacking of a VectorizedArray
template<typename Number , std::size_t width = internal::VectorizedArrayWidthSpecifier<Number>::max_width>
VectorizedArray< Number, width > make_vectorized_array (const Number &u)
 
template<typename VectorizedArrayType >
VectorizedArrayType make_vectorized_array (const typename VectorizedArrayType::value_type &u)
 
template<typename Number , std::size_t width>
void gather (VectorizedArray< Number, width > &out, const std::array< Number *, width > &ptrs, const unsigned int offset)
 
template<typename Number , std::size_t width>
void vectorized_load_and_transpose (const unsigned int n_entries, const Number *in, const unsigned int *offsets, VectorizedArray< Number, width > *out)
 
template<typename Number , std::size_t width>
void vectorized_load_and_transpose (const unsigned int n_entries, const std::array< Number *, width > &in, VectorizedArray< Number, width > *out)
 
template<typename Number , std::size_t width>
void vectorized_transpose_and_store (const bool add_into, const unsigned int n_entries, const VectorizedArray< Number, width > *in, const unsigned int *offsets, Number *out)
 
template<typename Number , std::size_t width>
void vectorized_transpose_and_store (const bool add_into, const unsigned int n_entries, const VectorizedArray< Number, width > *in, std::array< Number *, width > &out)
 
Arithmetic operations with VectorizedArray
template<typename Number , std::size_t width>
bool operator== (const VectorizedArray< Number, width > &lhs, const VectorizedArray< Number, width > &rhs)
 
template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator+ (const VectorizedArray< Number, width > &u, const VectorizedArray< Number, width > &v)
 
template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator- (const VectorizedArray< Number, width > &u, const VectorizedArray< Number, width > &v)
 
template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator* (const VectorizedArray< Number, width > &u, const VectorizedArray< Number, width > &v)
 
template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator/ (const VectorizedArray< Number, width > &u, const VectorizedArray< Number, width > &v)
 
template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator+ (const Number &u, const VectorizedArray< Number, width > &v)
 
template<std::size_t width>
VectorizedArray< float, width > operator+ (const double u, const VectorizedArray< float, width > &v)
 
template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator+ (const VectorizedArray< Number, width > &v, const Number &u)
 
template<std::size_t width>
VectorizedArray< float, width > operator+ (const VectorizedArray< float, width > &v, const double u)
 
template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator- (const Number &u, const VectorizedArray< Number, width > &v)
 
template<std::size_t width>
VectorizedArray< float, width > operator- (const double u, const VectorizedArray< float, width > &v)
 
template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator- (const VectorizedArray< Number, width > &v, const Number &u)
 
template<std::size_t width>
VectorizedArray< float, width > operator- (const VectorizedArray< float, width > &v, const double u)
 
template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator* (const Number &u, const VectorizedArray< Number, width > &v)
 
template<std::size_t width>
VectorizedArray< float, width > operator* (const double u, const VectorizedArray< float, width > &v)
 
template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator* (const VectorizedArray< Number, width > &v, const Number &u)
 
template<std::size_t width>
VectorizedArray< float, width > operator* (const VectorizedArray< float, width > &v, const double u)
 
template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator/ (const Number &u, const VectorizedArray< Number, width > &v)
 
template<std::size_t width>
VectorizedArray< float, width > operator/ (const double u, const VectorizedArray< float, width > &v)
 
template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator/ (const VectorizedArray< Number, width > &v, const Number &u)
 
template<std::size_t width>
VectorizedArray< float, width > operator/ (const VectorizedArray< float, width > &v, const double u)
 
template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator+ (const VectorizedArray< Number, width > &u)
 
template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator- (const VectorizedArray< Number, width > &u)
 
template<typename Number , std::size_t width>
std::ostream & operator<< (std::ostream &out, const VectorizedArray< Number, width > &p)
 

Ternary operations on VectorizedArray

enum class  SIMDComparison : int {
  equal , not_equal , less_than , less_than_or_equal ,
  greater_than , greater_than_or_equal
}
 
template<SIMDComparison predicate, typename Number >
Number compare_and_apply_mask (const Number &left, const Number &right, const Number &true_value, const Number &false_value)
 
template<SIMDComparison predicate, typename Number >
VectorizedArray< Number, 1 > compare_and_apply_mask (const VectorizedArray< Number, 1 > &left, const VectorizedArray< Number, 1 > &right, const VectorizedArray< Number, 1 > &true_value, const VectorizedArray< Number, 1 > &false_value)
 

Enumeration Type Documentation

◆ SIMDComparison

enum class SIMDComparison : int
strong

enum class encoding binary operations for a component-wise comparison of VectorizedArray data types.

Note
In case of SIMD vecorization (sse, avx, av512) we select the corresponding ordered, non-signalling (OQ) variants.
Enumerator
equal 
not_equal 
less_than 
less_than_or_equal 
greater_than 
greater_than_or_equal 

Definition at line 5354 of file vectorization.h.

Function Documentation

◆ make_vectorized_array() [1/2]

template<typename Number , std::size_t width = internal::VectorizedArrayWidthSpecifier<Number>::max_width>
VectorizedArray< Number, width > make_vectorized_array ( const Number &  u)
inline

Create a vectorized array that sets all entries in the array to the given scalar, i.e., broadcasts the scalar to all array elements.

Definition at line 764 of file vectorization.h.

◆ make_vectorized_array() [2/2]

template<typename VectorizedArrayType >
VectorizedArrayType make_vectorized_array ( const typename VectorizedArrayType::value_type &  u)
inline

Create a vectorized array of given type and broadcast the scalar value to all array elements.

Definition at line 780 of file vectorization.h.

◆ gather()

template<typename Number , std::size_t width>
void gather ( VectorizedArray< Number, width > &  out,
const std::array< Number *, width > &  ptrs,
const unsigned int  offset 
)
inline

Load size() data items from memory into the VectorizedArray out, starting at the given addresses and with given offset, each entry from the offset providing one element of the vectorized array.

This operation corresponds to the following code:

for (unsigned int v=0; v<VectorizedArray<Number>::size(); ++v)
out.data[v] = ptrs[v][offset];

Definition at line 807 of file vectorization.h.

◆ vectorized_load_and_transpose() [1/2]

template<typename Number , std::size_t width>
void vectorized_load_and_transpose ( const unsigned int  n_entries,
const Number *  in,
const unsigned int offsets,
VectorizedArray< Number, width > *  out 
)
inline

This method loads VectorizedArray::size() data streams from the given array in. The offsets to the input array are given by the array offsets. From each stream, n_entries are read. The data is then transposed and stored it into an array of VectorizedArray type. The output array out is expected to be an array of size n_entries. This method operates on plain arrays, so no checks for valid data access are made. It is the user's responsibility to ensure that the given arrays are valid according to the access layout below.

This operation corresponds to a transformation of an array-of-struct (input) into a struct-of-array (output) according to the following formula:

for (unsigned int i=0; i<n_entries; ++i)
for (unsigned int v=0; v<VectorizedArray<Number>::size(); ++v)
out[i][v] = in[offsets[v]+i];

A more optimized version of this code will be used for supported types.

This is the inverse operation to vectorized_transpose_and_store().

Definition at line 844 of file vectorization.h.

◆ vectorized_load_and_transpose() [2/2]

template<typename Number , std::size_t width>
void vectorized_load_and_transpose ( const unsigned int  n_entries,
const std::array< Number *, width > &  in,
VectorizedArray< Number, width > *  out 
)
inline

The same as above with the difference that an array of pointers are passed in as input argument in.

In analogy to the function above, one can consider that in+offset[v] is precomputed and passed as input argument.

However, this function can also be used if some function returns an array of pointers and no assumption can be made that they belong to the same array, i.e., they can have their origin in different memory allocations.

Definition at line 868 of file vectorization.h.

◆ vectorized_transpose_and_store() [1/2]

template<typename Number , std::size_t width>
void vectorized_transpose_and_store ( const bool  add_into,
const unsigned int  n_entries,
const VectorizedArray< Number, width > *  in,
const unsigned int offsets,
Number *  out 
)
inline

This method stores the vectorized arrays in transposed form into the given output array out with the given offsets offsets. This operation corresponds to a transformation of a struct-of-array (input) into an array- of-struct (output). This method operates on plain array, so no checks for valid data access are made. It is the user's responsibility to ensure that the given arrays are valid according to the access layout below.

This method assumes that the specified offsets do not overlap. Otherwise, the behavior is undefined in the vectorized case. It is the user's responsibility to make sure that the access does not overlap and avoid undefined behavior.

The argument add_into selects where the entries should only be written into the output arrays or the result should be added into the existing entries in the output. For add_into == false, the following code is assumed:

for (unsigned int i=0; i<n_entries; ++i)
for (unsigned int v=0; v<VectorizedArray<Number>::size(); ++v)
out[offsets[v]+i] = in[i][v];

For add_into == true, the code implements the following action:

for (unsigned int i=0; i<n_entries; ++i)
for (unsigned int v=0; v<VectorizedArray<Number>::size(); ++v)
out[offsets[v]+i] += in[i][v];

A more optimized version of this code will be used for supported types.

This is the inverse operation to vectorized_load_and_transpose().

Definition at line 919 of file vectorization.h.

◆ vectorized_transpose_and_store() [2/2]

template<typename Number , std::size_t width>
void vectorized_transpose_and_store ( const bool  add_into,
const unsigned int  n_entries,
const VectorizedArray< Number, width > *  in,
std::array< Number *, width > &  out 
)
inline

The same as above with the difference that an array of pointers are passed in as input argument out.

In analogy to the function above, one can consider that out+offset[v] is precomputed and passed as input argument.

However, this function can also be used if some function returns an array of pointers and no assumption can be made that they belong to the same array, i.e., they can have their origin in different memory allocations.

Definition at line 949 of file vectorization.h.

◆ operator==()

template<typename Number , std::size_t width>
bool operator== ( const VectorizedArray< Number, width > &  lhs,
const VectorizedArray< Number, width > &  rhs 
)
inline

Relational operator == for VectorizedArray

Definition at line 4994 of file vectorization.h.

◆ operator+() [1/6]

template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator+ ( const VectorizedArray< Number, width > &  u,
const VectorizedArray< Number, width > &  v 
)
inline

Addition of two vectorized arrays with operator +.

Definition at line 5012 of file vectorization.h.

◆ operator-() [1/6]

template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator- ( const VectorizedArray< Number, width > &  u,
const VectorizedArray< Number, width > &  v 
)
inline

Subtraction of two vectorized arrays with operator -.

Definition at line 5026 of file vectorization.h.

◆ operator*() [1/5]

template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator* ( const VectorizedArray< Number, width > &  u,
const VectorizedArray< Number, width > &  v 
)
inline

Multiplication of two vectorized arrays with operator *.

Definition at line 5040 of file vectorization.h.

◆ operator/() [1/5]

template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator/ ( const VectorizedArray< Number, width > &  u,
const VectorizedArray< Number, width > &  v 
)
inline

Division of two vectorized arrays with operator /.

Definition at line 5054 of file vectorization.h.

◆ operator+() [2/6]

template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator+ ( const Number &  u,
const VectorizedArray< Number, width > &  v 
)
inline

Addition of a scalar (expanded to a vectorized array with size() equal entries) and a vectorized array.

Definition at line 5069 of file vectorization.h.

◆ operator+() [3/6]

template<std::size_t width>
VectorizedArray< float, width > operator+ ( const double  u,
const VectorizedArray< float, width > &  v 
)
inline

Addition of a scalar (expanded to a vectorized array with size() equal entries) and a vectorized array in case the scalar is a double (needed in order to be able to write simple code with constants that are usually double numbers).

Definition at line 5085 of file vectorization.h.

◆ operator+() [4/6]

template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator+ ( const VectorizedArray< Number, width > &  v,
const Number &  u 
)
inline

Addition of a vectorized array and a scalar (expanded to a vectorized array with size() equal entries).

Definition at line 5099 of file vectorization.h.

◆ operator+() [5/6]

template<std::size_t width>
VectorizedArray< float, width > operator+ ( const VectorizedArray< float, width > &  v,
const double  u 
)
inline

Addition of a vectorized array and a scalar (expanded to a vectorized array with size() equal entries) in case the scalar is a double (needed in order to be able to write simple code with constants that are usually double numbers).

Definition at line 5114 of file vectorization.h.

◆ operator-() [2/6]

template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator- ( const Number &  u,
const VectorizedArray< Number, width > &  v 
)
inline

Subtraction of a vectorized array from a scalar (expanded to a vectorized array with size() equal entries).

Definition at line 5127 of file vectorization.h.

◆ operator-() [3/6]

template<std::size_t width>
VectorizedArray< float, width > operator- ( const double  u,
const VectorizedArray< float, width > &  v 
)
inline

Subtraction of a vectorized array from a scalar (expanded to a vectorized array with size() equal entries) in case the scalar is a double (needed in order to be able to write simple code with constants that are usually double numbers).

Definition at line 5143 of file vectorization.h.

◆ operator-() [4/6]

template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator- ( const VectorizedArray< Number, width > &  v,
const Number &  u 
)
inline

Subtraction of a scalar (expanded to a vectorized array with size() equal entries) from a vectorized array.

Definition at line 5157 of file vectorization.h.

◆ operator-() [5/6]

template<std::size_t width>
VectorizedArray< float, width > operator- ( const VectorizedArray< float, width > &  v,
const double  u 
)
inline

Subtraction of a scalar (expanded to a vectorized array with size() equal entries) from a vectorized array in case the scalar is a double (needed in order to be able to write simple code with constants that are usually double numbers).

Definition at line 5173 of file vectorization.h.

◆ operator*() [2/5]

template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator* ( const Number &  u,
const VectorizedArray< Number, width > &  v 
)
inline

Multiplication of a scalar (expanded to a vectorized array with size() equal entries) and a vectorized array.

Definition at line 5187 of file vectorization.h.

◆ operator*() [3/5]

template<std::size_t width>
VectorizedArray< float, width > operator* ( const double  u,
const VectorizedArray< float, width > &  v 
)
inline

Multiplication of a scalar (expanded to a vectorized array with size() equal entries) and a vectorized array in case the scalar is a double (needed in order to be able to write simple code with constants that are usually double numbers).

Definition at line 5203 of file vectorization.h.

◆ operator*() [4/5]

template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator* ( const VectorizedArray< Number, width > &  v,
const Number &  u 
)
inline

Multiplication of a vectorized array and a scalar (expanded to a vectorized array with size() equal entries).

Definition at line 5217 of file vectorization.h.

◆ operator*() [5/5]

template<std::size_t width>
VectorizedArray< float, width > operator* ( const VectorizedArray< float, width > &  v,
const double  u 
)
inline

Multiplication of a vectorized array and a scalar (expanded to a vectorized array with size() equal entries) in case the scalar is a double (needed in order to be able to write simple code with constants that are usually double numbers).

Definition at line 5232 of file vectorization.h.

◆ operator/() [2/5]

template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator/ ( const Number &  u,
const VectorizedArray< Number, width > &  v 
)
inline

Quotient between a scalar (expanded to a vectorized array with size() equal entries) and a vectorized array.

Definition at line 5245 of file vectorization.h.

◆ operator/() [3/5]

template<std::size_t width>
VectorizedArray< float, width > operator/ ( const double  u,
const VectorizedArray< float, width > &  v 
)
inline

Quotient between a scalar (expanded to a vectorized array with size() equal entries) and a vectorized array in case the scalar is a double (needed in order to be able to write simple code with constants that are usually double numbers).

Definition at line 5261 of file vectorization.h.

◆ operator/() [4/5]

template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator/ ( const VectorizedArray< Number, width > &  v,
const Number &  u 
)
inline

Quotient between a vectorized array and a scalar (expanded to a vectorized array with size() equal entries).

Definition at line 5275 of file vectorization.h.

◆ operator/() [5/5]

template<std::size_t width>
VectorizedArray< float, width > operator/ ( const VectorizedArray< float, width > &  v,
const double  u 
)
inline

Quotient between a vectorized array and a scalar (expanded to a vectorized array with size() equal entries) in case the scalar is a double (needed in order to be able to write simple code with constants that are usually double numbers).

Definition at line 5291 of file vectorization.h.

◆ operator+() [6/6]

template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator+ ( const VectorizedArray< Number, width > &  u)
inline

Unary operator + on a vectorized array.

Definition at line 5304 of file vectorization.h.

◆ operator-() [6/6]

template<typename Number , std::size_t width>
VectorizedArray< Number, width > operator- ( const VectorizedArray< Number, width > &  u)
inline

Unary operator - on a vectorized array.

Definition at line 5316 of file vectorization.h.

◆ operator<<()

template<typename Number , std::size_t width>
std::ostream & operator<< ( std::ostream &  out,
const VectorizedArray< Number, width > &  p 
)
inline

Output operator for vectorized array.

Definition at line 5329 of file vectorization.h.

◆ compare_and_apply_mask() [1/2]

template<SIMDComparison predicate, typename Number >
Number compare_and_apply_mask ( const Number &  left,
const Number &  right,
const Number &  true_value,
const Number &  false_value 
)
inline

Computes the vectorized equivalent of the following ternary operation:

(left OP right) ? true_value : false_value

where OP is a binary operator (such as =, !=, <, <=, >, and >=).

Such a computational idiom is useful as an alternative to branching whenever the control flow itself would depend on (computed) data. For example, in case of a scalar data type the statement (left < right) ? true_value : false_value could have been also implemented using an if-statement:

if (left < right)
result = true_value;
else
result = false_value;

This, however, is fundamentally impossible in case of vectorization because different decisions will be necessary on different vector entries (lanes) and the first variant (based on a ternary operator) has to be used instead:

result = compare_and_apply_mask<SIMDComparison::less_than>
(left, right, true_value, false_value);

Some more illustrative examples (that are less efficient than the dedicated std::max and std::abs overloads):

// std::max
const auto maximum = compare_and_apply_mask<SIMDComparison::greater_than>
(left, right, left, right);
// std::abs
const auto absolute = compare_and_apply_mask<SIMDComparison::less_than>
(left, VectorizedArray<double>(0.), -left, left);

More precisely, this function first computes a (boolean) mask that is the result of a binary operator OP applied to all elements of the VectorizedArray arguments left and right. The mask is then used to either select the corresponding component of true_value (if the binary operation equates to true), or false_value. The binary operator is encoded via the SIMDComparison template argument predicate.

In order to ease with generic programming approaches, the function provides overloads for all VectorizedArray<Number> variants as well as generic POD types such as double and float.

Note
For this function to work the binary operation has to be encoded via a SIMDComparison template argument predicate. Depending on it appropriate low-level machine instructions are generated replacing the call to compare_and_apply_mask. This also explains why predicate is a compile-time constant template parameter and not a constant function argument. In order to be able to emit the correct low-level instruction, the compiler has to know the comparison at compile time.

Definition at line 5439 of file vectorization.h.

◆ compare_and_apply_mask() [2/2]

template<SIMDComparison predicate, typename Number >
VectorizedArray< Number, 1 > compare_and_apply_mask ( const VectorizedArray< Number, 1 > &  left,
const VectorizedArray< Number, 1 > &  right,
const VectorizedArray< Number, 1 > &  true_value,
const VectorizedArray< Number, 1 > &  false_value 
)
inline

Specialization of above function for the non-vectorized VectorizedArray<Number, 1> variant.

Definition at line 5477 of file vectorization.h.