![]() |
Reference documentation for deal.II version GIT c00406db16 2023-09-28 18:00:02+00:00
|
#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.
Namespaces | |
internal | |
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) |
|
strong |
enum class encoding binary operations for a component-wise comparison of VectorizedArray data types.
OQ
) variants. Enumerator | |
---|---|
equal | |
not_equal | |
less_than | |
less_than_or_equal | |
greater_than | |
greater_than_or_equal |
Definition at line 5934 of file vectorization.h.
|
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 766 of file vectorization.h.
|
inline |
Create a vectorized array of given type and broadcast the scalar value to all array elements.
Definition at line 782 of file vectorization.h.
|
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:
Definition at line 809 of file vectorization.h.
|
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:
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 846 of file vectorization.h.
|
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 870 of file vectorization.h.
|
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 add_into == true
, the code implements the following action:
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 921 of file vectorization.h.
|
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 951 of file vectorization.h.
|
inline |
Relational operator == for VectorizedArray
Definition at line 5574 of file vectorization.h.
|
inline |
Addition of two vectorized arrays with operator +.
Definition at line 5592 of file vectorization.h.
|
inline |
Subtraction of two vectorized arrays with operator -.
Definition at line 5606 of file vectorization.h.
|
inline |
Multiplication of two vectorized arrays with operator *.
Definition at line 5620 of file vectorization.h.
|
inline |
Division of two vectorized arrays with operator /.
Definition at line 5634 of file vectorization.h.
|
inline |
Addition of a scalar (expanded to a vectorized array with size()
equal entries) and a vectorized array.
Definition at line 5649 of file vectorization.h.
|
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 5665 of file vectorization.h.
|
inline |
Addition of a vectorized array and a scalar (expanded to a vectorized array with size()
equal entries).
Definition at line 5679 of file vectorization.h.
|
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 5694 of file vectorization.h.
|
inline |
Subtraction of a vectorized array from a scalar (expanded to a vectorized array with size()
equal entries).
Definition at line 5707 of file vectorization.h.
|
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 5723 of file vectorization.h.
|
inline |
Subtraction of a scalar (expanded to a vectorized array with size()
equal entries) from a vectorized array.
Definition at line 5737 of file vectorization.h.
|
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 5753 of file vectorization.h.
|
inline |
Multiplication of a scalar (expanded to a vectorized array with size()
equal entries) and a vectorized array.
Definition at line 5767 of file vectorization.h.
|
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 5783 of file vectorization.h.
|
inline |
Multiplication of a vectorized array and a scalar (expanded to a vectorized array with size()
equal entries).
Definition at line 5797 of file vectorization.h.
|
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 5812 of file vectorization.h.
|
inline |
Quotient between a scalar (expanded to a vectorized array with size()
equal entries) and a vectorized array.
Definition at line 5825 of file vectorization.h.
|
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 5841 of file vectorization.h.
|
inline |
Quotient between a vectorized array and a scalar (expanded to a vectorized array with size()
equal entries).
Definition at line 5855 of file vectorization.h.
|
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 5871 of file vectorization.h.
|
inline |
Unary operator + on a vectorized array.
Definition at line 5884 of file vectorization.h.
|
inline |
Unary operator - on a vectorized array.
Definition at line 5896 of file vectorization.h.
|
inline |
Output operator for vectorized array.
Definition at line 5909 of file vectorization.h.
|
inline |
Computes the vectorized equivalent of the following ternary operation:
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:
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:
Some more illustrative examples (that are less efficient than the dedicated std::max
and std::abs
overloads):
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.
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 6019 of file vectorization.h.
|
inline |
Specialization of above function for the non-vectorized VectorizedArray<Number, 1> variant.
Definition at line 6057 of file vectorization.h.
inline ::VectorizedArray<Number, width> std::sin | ( | const ::VectorizedArray< Number, width > & | x | ) |
Compute the sine of a vectorized data field. The result is returned as vectorized array in the form {sin(x[0]), sin(x[1]), ..., sin(x[VectorizedArray::size()-1])}
.
Definition at line 6494 of file vectorization.h.
inline ::VectorizedArray<Number, width> std::cos | ( | const ::VectorizedArray< Number, width > & | x | ) |
Compute the cosine of a vectorized data field. The result is returned as vectorized array in the form {cos(x[0]), cos(x[1]), ..., cos(x[size()-1])}
.
Definition at line 6521 of file vectorization.h.
inline ::VectorizedArray<Number, width> std::tan | ( | const ::VectorizedArray< Number, width > & | x | ) |
Compute the tangent of a vectorized data field. The result is returned as vectorized array in the form {tan(x[0]), tan(x[1]), ..., tan(x[size()-1])}
.
Definition at line 6543 of file vectorization.h.
inline ::VectorizedArray<Number, width> std::exp | ( | const ::VectorizedArray< Number, width > & | x | ) |
Compute the exponential of a vectorized data field. The result is returned as vectorized array in the form {exp(x[0]), exp(x[1]), ..., exp(x[size()-1])}
.
Definition at line 6565 of file vectorization.h.
inline ::VectorizedArray<Number, width> std::log | ( | const ::VectorizedArray< Number, width > & | x | ) |
Compute the natural logarithm of a vectorized data field. The result is returned as vectorized array in the form {log(x[0]), log(x[1]), ..., log(x[size()-1])}
.
Definition at line 6587 of file vectorization.h.
inline ::VectorizedArray<Number, width> std::sqrt | ( | const ::VectorizedArray< Number, width > & | x | ) |
Compute the square root of a vectorized data field. The result is returned as vectorized array in the form {sqrt(x[0]), sqrt(x[1]), ..., sqrt(x[size()-1])}
.
Definition at line 6609 of file vectorization.h.
inline ::VectorizedArray<Number, width> std::pow | ( | const ::VectorizedArray< Number, width > & | x, |
const Number | p | ||
) |
Raises the given number x
to the power p
for a vectorized data field. The result is returned as vectorized array in the form {pow(x[0],p), pow(x[1],p), ..., pow(x[size()-1],p)}
.
Definition at line 6625 of file vectorization.h.
inline ::VectorizedArray<Number, width> std::pow | ( | const ::VectorizedArray< Number, width > & | x, |
const ::VectorizedArray< Number, width > & | p | ||
) |
Raises the given number x
to the power p
for a vectorized data field. The result is returned as vectorized array in the form {pow(x[0],p[0]), pow(x[1],p[1]), ..., pow(x[size()-1],p[size()-1])}
.
Definition at line 6648 of file vectorization.h.
inline ::VectorizedArray<Number, width> std::abs | ( | const ::VectorizedArray< Number, width > & | x | ) |
Compute the absolute value (modulus) of a vectorized data field. The result is returned as vectorized array in the form {abs(x[0]), abs(x[1]), ..., abs(x[size()-1])}
.
Definition at line 6671 of file vectorization.h.
inline ::VectorizedArray<Number, width> std::max | ( | const ::VectorizedArray< Number, width > & | x, |
const ::VectorizedArray< Number, width > & | y | ||
) |
Compute the componentwise maximum of two vectorized data fields. The result is returned as vectorized array in the form {max(x[0],y[0]), max(x[1],y[1]), ...}
.
Definition at line 6687 of file vectorization.h.
inline ::VectorizedArray<Number, width> std::min | ( | const ::VectorizedArray< Number, width > & | x, |
const ::VectorizedArray< Number, width > & | y | ||
) |
Compute the componentwise minimum of two vectorized data fields. The result is returned as vectorized array in the form {min(x[0],y[0]), min(x[1],y[1]), ...}
.
Definition at line 6704 of file vectorization.h.