// for safety, also check that __AVX512F__ is defined in case the user manually
// set some conflicting compile flags which prevent compilation
-#if DEAL_II_COMPILER_VECTORIZATION_LEVEL == 3 && defined(__AVX512F__)
+#if DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 3 && defined(__AVX512F__)
/**
* Specialization of VectorizedArray class for double and AVX-512.
VectorizedArray &
operator = (const double x)
{
- data = _mm512_set_pd(x, x, x, x, x, x, x, x);
+ data = _mm512_set1_pd(x);
return *this;
}
operator += (const VectorizedArray &vec)
{
// if the compiler supports vector arithmetics, we can simply use +=
- // operator on the given data type. Otherwise, we need to use the built-in
+ // operator on the given data type. this allows the compiler to combine
+ // additions with multiplication (fused multiply-add) if those
+ // instructions are available. Otherwise, we need to use the built-in
// intrinsic command for __m512d
#ifdef DEAL_II_COMPILER_USE_VECTOR_ARITHMETICS
data += vec.data;
{
// to compute the absolute value, perform bitwise andnot with -0. This
// will leave all value and exponent bits unchanged but force the sign
- // value to +.
- __m256d mask = _mm512_set_pd (-0., -0., -0., -0., -0., -0., -0., -0.);
+ // value to +. As opposed to AVX and SSE2, there is no andnot operation on
+ // double data types in AVX-512. Thus, need to separately work with
+ // 256-bit data fields.
+ __m256d mask = _mm256_set1_pd (-0.);
+ __m256d in[2];
+ in[0] = *(reinterpret_cast<__m256d *>(&data));
+ in[0] = _m256_andnot_pd(mask, in[0]);
+ in[1] = *(reinterpret_cast<__m256d *>(&data)+1);
+ in[1] = _m256_andnot_pd(mask, in[1]);
VectorizedArray res;
- res.data = _mm512_andnot_pd(mask, data);
+ res.data = _mm512_loadu_pd(reinterpret_cast<double *>(&in[0]));
return res;
}
VectorizedArray &
operator = (const float x)
{
- data = _mm512_set_ps(x, x, x, x, x, x, x, x, x, x, x, x, x, x, x, x);
+ data = _mm512_set1_ps(x);
return *this;
}
VectorizedArray &
operator += (const VectorizedArray &vec)
{
+ // if the compiler supports vector arithmetics, we can simply use +=
+ // operator on the given data type. this allows the compiler to combine
+ // additions with multiplication (fused multiply-add) if those
+ // instructions are available. Otherwise, we need to use the built-in
+ // intrinsic command for __m512d
#ifdef DEAL_II_COMPILER_USE_VECTOR_ARITHMETICS
data += vec.data;
#else
VectorizedArray
get_abs () const
{
- // to compute the absolute value, perform
- // bitwise andnot with -0. This will leave all
- // value and exponent bits unchanged but force
- // the sign value to +.
- __m256 mask = _mm512_setzero_ps (-0.f, -0.f, -0.f, -0.f, -0.f, -0.f, -0.f, -0.f,
- -0.f, -0.f, -0.f, -0.f, -0.f, -0.f, -0.f, -0.f);
+ // to compute the absolute value, perform bitwise andnot with -0. This
+ // will leave all value and exponent bits unchanged but force the sign
+ // value to +. As opposed to AVX and SSE, there is no andnot operation on
+ // double data types in AVX-512. Thus, need to separately work with
+ // 256-bit data fields.
+ __m256 mask = _mm256_set1_ps (-0.f);
+ __m256 in[2];
+ in[0] = *(reinterpret_cast<__m256 *>(&data));
+ in[0] = _m256_andnot_ps(mask, in[0]);
+ in[1] = *(reinterpret_cast<__m256 *>(&data)+1);
+ in[1] = _m256_andnot_ps(mask, in[1]);
VectorizedArray res;
- res.data = _mm512_andnot_ps(mask, data);
+ res.data = _mm512_loadu_ps(reinterpret_cast<float *>(&in[0]));
return res;
}
-#elif DEAL_II_COMPILER_VECTORIZATION_LEVEL == 2 && defined(__AVX__)
+#elif DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 2 && defined(__AVX__)
/**
* Specialization of VectorizedArray class for double and AVX.
{
public:
/**
- * This gives the number of vectors collected
- * in this class.
+ * This gives the number of vectors collected in this class.
*/
static const unsigned int n_array_elements = 4;
/**
- * This function can be used to set all data
- * fields to a given scalar.
+ * This function can be used to set all data fields to a given scalar.
*/
VectorizedArray &
operator = (const double x)
VectorizedArray &
operator += (const VectorizedArray &vec)
{
- // if the compiler supports vector
- // arithmetics, we can simply use += operator
- // on the given data type. Otherwise, we need
- // to use the built-in intrinsic command for
- // __m256d
+ // if the compiler supports vector arithmetics, we can simply use +=
+ // operator on the given data type. this allows the compiler to combine
+ // additions with multiplication (fused multiply-add) if those
+ // instructions are available. Otherwise, we need to use the built-in
+ // intrinsic command for __m256d
#ifdef DEAL_II_COMPILER_USE_VECTOR_ARITHMETICS
data += vec.data;
#else
}
/**
- * Actual data field. Since this class
- * represents a POD data type, it remains
- * public.
+ * Actual data field. Since this class represents a POD data type, it
+ * remains public.
*/
__m256d data;
private:
/**
- * Returns the square root of this field. Not
- * for use in user code. Use sqrt(x) instead.
+ * Returns the square root of this field. Not for use in user code. Use
+ * sqrt(x) instead.
*/
VectorizedArray
get_sqrt () const
}
/**
- * Returns the absolute value of this
- * field. Not for use in user code. Use
+ * Returns the absolute value of this field. Not for use in user code. Use
* abs(x) instead.
*/
VectorizedArray
get_abs () const
{
- // to compute the absolute value, perform
- // bitwise andnot with -0. This will leave all
- // value and exponent bits unchanged but force
- // the sign value to +.
+ // to compute the absolute value, perform bitwise andnot with -0. This
+ // will leave all value and exponent bits unchanged but force the sign
+ // value to +.
__m256d mask = _mm256_set1_pd (-0.);
VectorizedArray res;
res.data = _mm256_andnot_pd(mask, data);
}
/**
- * Returns the component-wise maximum of this
- * field and another one. Not for use in user
- * code. Use max(x,y) instead.
+ * Returns the component-wise maximum of this field and another one. Not for
+ * use in user code. Use max(x,y) instead.
*/
VectorizedArray
get_max (const VectorizedArray &other) const
}
/**
- * Returns the component-wise minimum of this
- * field and another one. Not for use in user
- * code. Use min(x,y) instead.
+ * Returns the component-wise minimum of this field and another one. Not for
+ * use in user code. Use min(x,y) instead.
*/
VectorizedArray
get_min (const VectorizedArray &other) const
{
public:
/**
- * This gives the number of vectors collected
- * in this class.
+ * This gives the number of vectors collected in this class.
*/
static const unsigned int n_array_elements = 8;
/**
- * This function can be used to set all data
- * fields to a given scalar.
+ * This function can be used to set all data fields to a given scalar.
*/
VectorizedArray &
operator = (const float x)
VectorizedArray &
operator += (const VectorizedArray &vec)
{
+ // if the compiler supports vector arithmetics, we can simply use +=
+ // operator on the given data type. this allows the compiler to combine
+ // additions with multiplication (fused multiply-add) if those
+ // instructions are available. Otherwise, we need to use the built-in
+ // intrinsic command for __m256d
#ifdef DEAL_II_COMPILER_USE_VECTOR_ARITHMETICS
data += vec.data;
#else
}
/**
- * Actual data field. Since this class
- * represents a POD data type, it remains
- * public.
+ * Actual data field. Since this class represents a POD data type, it
+ * remains public.
*/
__m256 data;
private:
/**
- * Returns the square root of this field. Not
- * for use in user code. Use sqrt(x) instead.
+ * Returns the square root of this field. Not for use in user code. Use
+ * sqrt(x) instead.
*/
VectorizedArray
get_sqrt () const
res.data = _mm256_sqrt_ps(data);
return res;
}
+
/**
- * Returns the absolute value of this
- * field. Not for use in user code. Use
+ * Returns the absolute value of this field. Not for use in user code. Use
* abs(x) instead.
*/
VectorizedArray
get_abs () const
{
- // to compute the absolute value, perform
- // bitwise andnot with -0. This will leave all
- // value and exponent bits unchanged but force
- // the sign value to +.
+ // to compute the absolute value, perform bitwise andnot with -0. This
+ // will leave all value and exponent bits unchanged but force the sign
+ // value to +.
__m256 mask = _mm256_set1_ps (-0.f);
VectorizedArray res;
res.data = _mm256_andnot_ps(mask, data);
}
/**
- * Returns the component-wise maximum of this
- * field and another one. Not for use in user
- * code. Use max(x,y) instead.
+ * Returns the component-wise maximum of this field and another one. Not for
+ * use in user code. Use max(x,y) instead.
*/
VectorizedArray
get_max (const VectorizedArray &other) const
}
/**
- * Returns the component-wise minimum of this
- * field and another one. Not for use in user
- * code. Use min(x,y) instead.
+ * Returns the component-wise minimum of this field and another one. Not for
+ * use in user code. Use min(x,y) instead.
*/
VectorizedArray
get_min (const VectorizedArray &other) const
{
public:
/**
- * This gives the number of vectors collected
- * in this class.
+ * This gives the number of vectors collected in this class.
*/
static const unsigned int n_array_elements = 2;
/**
- * This function can be used to set all data
- * fields to a given scalar.
+ * This function can be used to set all data fields to a given scalar.
*/
-
VectorizedArray &
operator = (const double x)
{
}
/**
- * Actual data field. Since this class
- * represents a POD data type, it remains
- * public.
+ * Actual data field. Since this class represents a POD data type, it
+ * remains public.
*/
__m128d data;
private:
/**
- * Returns the square root of this field. Not
- * for use in user code. Use sqrt(x) instead.
+ * Returns the square root of this field. Not for use in user code. Use
+ * sqrt(x) instead.
*/
VectorizedArray
get_sqrt () const
}
/**
- * Returns the absolute value of this
- * field. Not for use in user code. Use abs(x)
- * instead.
+ * Returns the absolute value of this field. Not for use in user code. Use
+ * abs(x) instead.
*/
VectorizedArray
get_abs () const
}
/**
- * Returns the component-wise maximum of this
- * field and another one. Not for use in user
- * code. Use max(x,y) instead.
+ * Returns the component-wise maximum of this field and another one. Not for
+ * use in user code. Use max(x,y) instead.
*/
VectorizedArray
get_max (const VectorizedArray &other) const
}
/**
- * Returns the component-wise minimum of this
- * field and another one. Not for use in user
- * code. Use min(x,y) instead.
+ * Returns the component-wise minimum of this field and another one. Not for
+ * use in user code. Use min(x,y) instead.
*/
VectorizedArray
get_min (const VectorizedArray &other) const
{
public:
/**
- * This gives the number of vectors collected
- * in this class.
+ * This gives the number of vectors collected in this class.
*/
static const unsigned int n_array_elements = 4;
/**
- * This function can be used to set all data
- * fields to a given scalar.
+ * This function can be used to set all data fields to a given scalar.
*/
VectorizedArray &
}
/**
- * Actual data field. Since this class
- * represents a POD data type, it remains
- * public.
+ * Actual data field. Since this class represents a POD data type, it
+ * remains public.
*/
__m128 data;
private:
/**
- * Returns the square root of this field. Not
- * for use in user code. Use sqrt(x) instead.
+ * Returns the square root of this field. Not for use in user code. Use
+ * sqrt(x) instead.
*/
VectorizedArray
get_sqrt () const
}
/**
- * Returns the absolute value of this
- * field. Not for use in user code. Use
+ * Returns the absolute value of this field. Not for use in user code. Use
* abs(x) instead.
*/
VectorizedArray
get_abs () const
{
- // to compute the absolute value, perform
- // bitwise andnot with -0. This will leave all
- // value and exponent bits unchanged but force
- // the sign value to +.
+ // to compute the absolute value, perform bitwise andnot with -0. This
+ // will leave all value and exponent bits unchanged but force the sign
+ // value to +.
__m128 mask = _mm_set1_ps (-0.f);
VectorizedArray res;
res.data = _mm_andnot_ps(mask, data);
}
/**
- * Returns the component-wise maximum of this
- * field and another one. Not for use in user
- * code. Use max(x,y) instead.
+ * Returns the component-wise maximum of this field and another one. Not for
+ * use in user code. Use max(x,y) instead.
*/
VectorizedArray
get_max (const VectorizedArray &other) const
}
/**
- * Returns the component-wise minimum of this
- * field and another one. Not for use in user
- * code. Use min(x,y) instead.
+ * Returns the component-wise minimum of this field and another one. Not for
+ * use in user code. Use min(x,y) instead.
*/
VectorizedArray
get_min (const VectorizedArray &other) const
* when compiling deal.II on 64-bit operating systems. On Intel Sandy Bridge
* processors and newer or AMD Bulldozer processors and newer, four doubles
* and eight floats are used when deal.II is configured e.g. using gcc with
- * --with-cpu=native or --with-cpu=corei7-avx.
+ * --with-cpu=native or --with-cpu=corei7-avx. On compilations with AVX-512
+ * support, eight doubles and sixteen floats are used.
*
* This behavior of this class is made similar to the basic data types
* double and float. The definition of a vectorized array does not
* VectorizedArray<double>()</tt>, the compiler will take care of data
* alignment automatically. However, when allocating a long vector of
* VectorizedArray<double> data, one needs to respect these rules. Use the
- * class AlignedVector for this purpose. It is a class very similar to
+ * class AlignedVector or data containers based on AlignedVector (such as
+ * Table) for this purpose. It is a class very similar to
* std::vector otherwise but always makes sure that data is correctly
* aligned.
*
{
public:
/**
- * This gives the number of vectors collected
- * in this class.
+ * This gives the number of vectors collected in this class.
*/
static const unsigned int n_array_elements = 1;
- // POD means that there should be no
- // user-defined constructors, destructors and
- // copy functions (the standard is somewhat
- // relaxed in C++2011, though).
+ // POD means that there should be no user-defined constructors, destructors
+ // and copy functions (the standard is somewhat relaxed in C++2011, though).
/**
- * This function assigns a scalar to this
- * class.
+ * This function assigns a scalar to this class.
*/
-
VectorizedArray &
operator = (const Number scalar)
{
}
/**
- * Access operator (only valid with component
- * 0)
+ * Access operator (only valid with component 0)
*/
Number &
operator [] (const unsigned int comp)
}
/**
- * Constant access operator (only valid with
- * component 0)
+ * Constant access operator (only valid with component 0)
*/
const Number &
operator [] (const unsigned int comp) const
}
/**
- * Actual data field. Since this class
- * represents a POD data type, it is declared
- * public.
+ * Actual data field. Since this class represents a POD data type, it is
+ * declared public.
*/
Number data;
private:
/**
- * Returns the square root of this field. Not
- * for use in user code. Use sqrt(x) instead.
+ * Returns the square root of this field. Not for use in user code. Use
+ * sqrt(x) instead.
*/
VectorizedArray
get_sqrt () const
}
/**
- * Returns the absolute value of this
- * field. Not for use in user code. Use
+ * Returns the absolute value of this field. Not for use in user code. Use
* abs(x) instead.
*/
VectorizedArray
}
/**
- * Returns the component-wise maximum of this
- * field and another one. Not for use in user
- * code. Use max(x,y) instead.
+ * Returns the component-wise maximum of this field and another one. Not for
+ * use in user code. Use max(x,y) instead.
*/
VectorizedArray
get_max (const VectorizedArray &other) const
}
/**
- * Returns the component-wise minimum of this
- * field and another one. Not for use in user
- * code. Use min(x,y) instead.
+ * Returns the component-wise minimum of this field and another one. Not for
+ * use in user code. Use min(x,y) instead.
*/
VectorizedArray
get_min (const VectorizedArray &other) const
sin (const ::dealii::VectorizedArray<Number> &x)
{
// put values in an array and later read in that array with an unaligned
- // read. This should safe some instructions as compared to directly
+ // read. This should save some instructions as compared to directly
// setting the individual elements and also circumvents a compiler
- // optimization bug in gcc-4.6 (see also deal.II developers list from
- // April 2014, topic "matrix_free/step-48 Test").
+ // optimization bug in gcc-4.6 with SSE2 (see also deal.II developers list
+ // from April 2014, topic "matrix_free/step-48 Test").
Number values[::dealii::VectorizedArray<Number>::n_array_elements];
for (unsigned int i=0; i<dealii::VectorizedArray<Number>::n_array_elements; ++i)
values[i] = std::sin(x[i]);