From: peterrum Date: Mon, 1 Jul 2019 13:09:52 +0000 (+0200) Subject: Extend VectorizedArray with a second template argument X-Git-Tag: v9.2.0-rc1~1398^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=f6ec327adb5ecb50fbf2998f1313566f07215ea6;p=dealii.git Extend VectorizedArray with a second template argument --- diff --git a/doc/news/changes/minor/20190704PeterMunch b/doc/news/changes/minor/20190704PeterMunch new file mode 100644 index 0000000000..2d5f56b552 --- /dev/null +++ b/doc/news/changes/minor/20190704PeterMunch @@ -0,0 +1,6 @@ +New: A second template parameter has been added to VectorizedArray. This template +parameter optionally controls the number of vectorization lanes to be used, i.e. +which instruction set architecture extension to be used, at the time the application +is built. +
+(Peter Munch, 2019/07/04) diff --git a/include/deal.II/base/memory_consumption.h b/include/deal.II/base/memory_consumption.h index 4d3fcd4d9d..85abcd1f35 100644 --- a/include/deal.II/base/memory_consumption.h +++ b/include/deal.II/base/memory_consumption.h @@ -30,12 +30,6 @@ DEAL_II_NAMESPACE_OPEN - -// forward declaration -template -class VectorizedArray; - - /** * This namespace provides functions helping to determine the amount of memory * used by objects. The goal is not necessarily to give the amount of memory @@ -134,9 +128,9 @@ namespace MemoryConsumption * Determine the amount of memory in bytes consumed by a * VectorizedArray variable. */ - template + template inline std::size_t - memory_consumption(const VectorizedArray &); + memory_consumption(const VectorizedArray &); /** * Determine an estimate of the amount of memory in bytes consumed by a @@ -302,11 +296,11 @@ namespace MemoryConsumption - template + template inline std::size_t - memory_consumption(const VectorizedArray &) + memory_consumption(const VectorizedArray &) { - return sizeof(VectorizedArray); + return sizeof(VectorizedArray); } diff --git a/include/deal.II/base/numbers.h b/include/deal.II/base/numbers.h index 157c5ec8be..a2e1b1de74 100644 --- a/include/deal.II/base/numbers.h +++ b/include/deal.II/base/numbers.h @@ -37,8 +37,55 @@ DEAL_II_NAMESPACE_OPEN +namespace internal +{ + template + struct VectorizedArrayWidthSpecifier + { + static const unsigned int max_width = 1; + }; + + template <> + struct VectorizedArrayWidthSpecifier + { + static const unsigned int max_width = +#if DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 1 && defined(__ALTIVEC__) + 2; +#elif DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 3 && defined(__AVX512F__) + 8; +#elif DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 2 && defined(__AVX__) + 4; +#elif DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 1 && defined(__SSE2__) + 2; +#else + 1; +#endif + }; + + template <> + struct VectorizedArrayWidthSpecifier + { + static const unsigned int max_width = +#if DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 1 && defined(__ALTIVEC__) + 4; +#elif DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 3 && defined(__AVX512F__) + 16; +#elif DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 2 && defined(__AVX__) + 8; +#elif DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 1 && defined(__SSE2__) + 4; +#else + 1; +#endif + }; + + +} // namespace internal + // forward declarations to support abs or sqrt operations on VectorizedArray -template +template ::max_width> class VectorizedArray; template struct EnableIfScalar; @@ -62,38 +109,38 @@ DEAL_II_NAMESPACE_CLOSE namespace std { - template - DEAL_II_ALWAYS_INLINE ::dealii::VectorizedArray - sqrt(const ::dealii::VectorizedArray &); - template - DEAL_II_ALWAYS_INLINE ::dealii::VectorizedArray - abs(const ::dealii::VectorizedArray &); - template - DEAL_II_ALWAYS_INLINE ::dealii::VectorizedArray - max(const ::dealii::VectorizedArray &, - const ::dealii::VectorizedArray &); - template - DEAL_II_ALWAYS_INLINE ::dealii::VectorizedArray - min(const ::dealii::VectorizedArray &, - const ::dealii::VectorizedArray &); - template - ::dealii::VectorizedArray - pow(const ::dealii::VectorizedArray &, const Number p); - template - ::dealii::VectorizedArray - sin(const ::dealii::VectorizedArray &); - template - ::dealii::VectorizedArray - cos(const ::dealii::VectorizedArray &); - template - ::dealii::VectorizedArray - tan(const ::dealii::VectorizedArray &); - template - ::dealii::VectorizedArray - exp(const ::dealii::VectorizedArray &); - template - ::dealii::VectorizedArray - log(const ::dealii::VectorizedArray &); + template + DEAL_II_ALWAYS_INLINE ::dealii::VectorizedArray + sqrt(const ::dealii::VectorizedArray &); + template + DEAL_II_ALWAYS_INLINE ::dealii::VectorizedArray + abs(const ::dealii::VectorizedArray &); + template + DEAL_II_ALWAYS_INLINE ::dealii::VectorizedArray + max(const ::dealii::VectorizedArray &, + const ::dealii::VectorizedArray &); + template + DEAL_II_ALWAYS_INLINE ::dealii::VectorizedArray + min(const ::dealii::VectorizedArray &, + const ::dealii::VectorizedArray &); + template + ::dealii::VectorizedArray + pow(const ::dealii::VectorizedArray &, const Number p); + template + ::dealii::VectorizedArray + sin(const ::dealii::VectorizedArray &); + template + ::dealii::VectorizedArray + cos(const ::dealii::VectorizedArray &); + template + ::dealii::VectorizedArray + tan(const ::dealii::VectorizedArray &); + template + ::dealii::VectorizedArray + exp(const ::dealii::VectorizedArray &); + template + ::dealii::VectorizedArray + log(const ::dealii::VectorizedArray &); } // namespace std DEAL_II_NAMESPACE_OPEN diff --git a/include/deal.II/base/tensor.h b/include/deal.II/base/tensor.h index d9f8f4c676..dec9a56867 100644 --- a/include/deal.II/base/tensor.h +++ b/include/deal.II/base/tensor.h @@ -45,8 +45,6 @@ template class Tensor; template class Vector; -template -class VectorizedArray; #ifndef DOXYGEN // Overload invalid tensor types of negative rank that come up during @@ -740,27 +738,27 @@ namespace internal } }; - template - struct NumberType>> + template + struct NumberType>> { - static constexpr const Tensor> & - value(const Tensor> &t) + static constexpr const Tensor> & + value(const Tensor> &t) { return t; } - static DEAL_II_CONSTEXPR Tensor> + static DEAL_II_CONSTEXPR Tensor> value(const T &t) { - Tensor> tmp; - tmp = internal::NumberType>::value(t); + Tensor> tmp; + tmp = internal::NumberType>::value(t); return tmp; } - static DEAL_II_CONSTEXPR Tensor> - value(const VectorizedArray &t) + static DEAL_II_CONSTEXPR Tensor> + value(const VectorizedArray &t) { - Tensor> tmp; + Tensor> tmp; tmp = t; return tmp; } diff --git a/include/deal.II/base/vectorization.h b/include/deal.II/base/vectorization.h index 19f9182a76..0ab3c74eb4 100644 --- a/include/deal.II/base/vectorization.h +++ b/include/deal.II/base/vectorization.h @@ -88,19 +88,19 @@ namespace internal * https://github.com/dealii/dealii/pull/3967 . Also see numbers.h for other * specializations. */ - template - struct NumberType> + template + struct NumberType> { - static const VectorizedArray & - value(const VectorizedArray &t) + static const VectorizedArray & + value(const VectorizedArray &t) { return t; } - static VectorizedArray + static VectorizedArray value(const T &t) { - VectorizedArray tmp; + VectorizedArray tmp; tmp = t; return tmp; } @@ -111,8 +111,8 @@ namespace internal // Enable the EnableIfScalar type trait for VectorizedArray such // that it can be used as a Number type in Tensor, etc. -template -struct EnableIfScalar> +template +struct EnableIfScalar> { using type = VectorizedArray::type>; }; @@ -167,9 +167,42 @@ struct EnableIfScalar> * similar to std::vector otherwise but always makes sure that data is * correctly aligned. * - * @author Katharina Kormann, Martin Kronbichler, 2010, 2011 + * The user can explicitly control the width of a particular instruction set + * architecture (ISA) extension by specifying the number of lanes via the second + * template parameter of this wrapper class. For example on Intel Skylake, + * you have the following options for the data type double: + * - VectorizedArray // no vectorization (auto-optimization) + * - VectorizedArray // SSE2 + * - VectorizedArray // AVX + * - VectorizedArray // AVX-512 (default) + * + * and for Intel Sandy Bridge, Haswell, Broadwell, AMD Bulldozer and Zen/Ryzen: + * - VectorizedArray // no vectorization (auto-optimization) + * - VectorizedArray // SSE2 + * - VectorizedArray // AVX (default) + * + * and for processors with AltiVec support: + * - VectorizedArray + * - VectorizedArray + * + * for older x86 processors or in case no processor-specific compilation flags + * were added (i.e., without `-D CMAKE_CXX_FLAGS=-march=native` or similar + * flags): + * - VectorizedArray // no vectorization (auto-optimization) + * - VectorizedArray // SSE2 + * + * Similar considerations also apply to the data type `float`. + * + * Wrongly selecting the width i.e width=3 or width>8 for Skylake leads to + * a static assert. + * + * @tparam Number underlying data type + * @tparam width vector length (optional; if not set, the maximal width of the + * architecture is used) + * + * @author Katharina Kormann, Martin Kronbichler, Peter Munch, 2010, 2011, 2019 */ -template +template class VectorizedArray { public: @@ -180,6 +213,9 @@ public: */ static const unsigned int n_array_elements = 1; + static_assert(width == n_array_elements, + "You specified an illegal width that is not supported."); + // POD means that there should be no user-defined constructors, destructors // and copy functions (the standard is somewhat relaxed in C++2011, though). @@ -221,7 +257,7 @@ public: */ DEAL_II_ALWAYS_INLINE VectorizedArray & - operator+=(const VectorizedArray &vec) + operator+=(const VectorizedArray &vec) { data += vec.data; return *this; @@ -232,7 +268,7 @@ public: */ DEAL_II_ALWAYS_INLINE VectorizedArray & - operator-=(const VectorizedArray &vec) + operator-=(const VectorizedArray &vec) { data -= vec.data; return *this; @@ -243,7 +279,7 @@ public: */ DEAL_II_ALWAYS_INLINE VectorizedArray & - operator*=(const VectorizedArray &vec) + operator*=(const VectorizedArray &vec) { data *= vec.data; return *this; @@ -254,7 +290,7 @@ public: */ DEAL_II_ALWAYS_INLINE VectorizedArray & - operator/=(const VectorizedArray &vec) + operator/=(const VectorizedArray &vec) { data /= vec.data; return *this; @@ -437,23 +473,27 @@ private: /** * Make a few functions friends. */ - template - friend VectorizedArray - std::sqrt(const VectorizedArray &); - template - friend VectorizedArray - std::abs(const VectorizedArray &); - template - friend VectorizedArray - std::max(const VectorizedArray &, const VectorizedArray &); - template - friend VectorizedArray - std::min(const VectorizedArray &, const VectorizedArray &); + template + friend VectorizedArray + std::sqrt(const VectorizedArray &); + template + friend VectorizedArray + std::abs(const VectorizedArray &); + template + friend VectorizedArray + std::max(const VectorizedArray &, + const VectorizedArray &); + template + friend VectorizedArray + std::min(const VectorizedArray &, + const VectorizedArray &); }; + + // We need to have a separate declaration for static const members -template -const unsigned int VectorizedArray::n_array_elements; +template +const unsigned int VectorizedArray::n_array_elements; @@ -463,11 +503,13 @@ const unsigned int VectorizedArray::n_array_elements; * * @relatesalso VectorizedArray */ -template -inline DEAL_II_ALWAYS_INLINE VectorizedArray +template < + typename Number, + int width = internal::VectorizedArrayWidthSpecifier::max_width> +inline DEAL_II_ALWAYS_INLINE VectorizedArray make_vectorized_array(const Number &u) { - VectorizedArray result; + VectorizedArray result; result = u; return result; } @@ -499,15 +541,17 @@ inline DEAL_II_ALWAYS_INLINE VectorizedArray * * @relatesalso VectorizedArray */ -template +template inline DEAL_II_ALWAYS_INLINE void -vectorized_load_and_transpose(const unsigned int n_entries, - const Number * in, - const unsigned int * offsets, - VectorizedArray *out) +vectorized_load_and_transpose(const unsigned int n_entries, + const Number * in, + const unsigned int * offsets, + VectorizedArray *out) { for (unsigned int i = 0; i < n_entries; ++i) - for (unsigned int v = 0; v < VectorizedArray::n_array_elements; ++v) + for (unsigned int v = 0; + v < VectorizedArray::n_array_elements; + ++v) out[i][v] = in[offsets[v] + i]; } @@ -551,28 +595,29 @@ vectorized_load_and_transpose(const unsigned int n_entries, * * @relatesalso VectorizedArray */ -template +template inline DEAL_II_ALWAYS_INLINE void -vectorized_transpose_and_store(const bool add_into, - const unsigned int n_entries, - const VectorizedArray *in, - const unsigned int * offsets, - Number * out) +vectorized_transpose_and_store(const bool add_into, + const unsigned int n_entries, + const VectorizedArray *in, + const unsigned int * offsets, + Number * out) { if (add_into) for (unsigned int i = 0; i < n_entries; ++i) - for (unsigned int v = 0; v < VectorizedArray::n_array_elements; + for (unsigned int v = 0; + v < VectorizedArray::n_array_elements; ++v) out[offsets[v] + i] += in[i][v]; else for (unsigned int i = 0; i < n_entries; ++i) - for (unsigned int v = 0; v < VectorizedArray::n_array_elements; + for (unsigned int v = 0; + v < VectorizedArray::n_array_elements; ++v) out[offsets[v] + i] = in[i][v]; } - // for safety, also check that __AVX512F__ is defined in case the user manually // set some conflicting compile flags which prevent compilation @@ -582,7 +627,7 @@ vectorized_transpose_and_store(const bool add_into, * Specialization of VectorizedArray class for double and AVX-512. */ template <> -class VectorizedArray +class VectorizedArray { public: /** @@ -848,18 +893,20 @@ private: /** * Make a few functions friends. */ - template - friend VectorizedArray - std::sqrt(const VectorizedArray &); - template - friend VectorizedArray - std::abs(const VectorizedArray &); - template - friend VectorizedArray - std::max(const VectorizedArray &, const VectorizedArray &); - template - friend VectorizedArray - std::min(const VectorizedArray &, const VectorizedArray &); + template + friend VectorizedArray + std::sqrt(const VectorizedArray &); + template + friend VectorizedArray + std::abs(const VectorizedArray &); + template + friend VectorizedArray + std::max(const VectorizedArray &, + const VectorizedArray &); + template + friend VectorizedArray + std::min(const VectorizedArray &, + const VectorizedArray &); }; @@ -869,10 +916,10 @@ private: */ template <> inline DEAL_II_ALWAYS_INLINE void -vectorized_load_and_transpose(const unsigned int n_entries, - const double * in, - const unsigned int * offsets, - VectorizedArray *out) +vectorized_load_and_transpose(const unsigned int n_entries, + const double * in, + const unsigned int * offsets, + VectorizedArray *out) { // do not do full transpose because the code is long and will most // likely not pay off because many processors have two load units @@ -914,11 +961,11 @@ vectorized_load_and_transpose(const unsigned int n_entries, */ template <> inline DEAL_II_ALWAYS_INLINE void -vectorized_transpose_and_store(const bool add_into, - const unsigned int n_entries, - const VectorizedArray *in, - const unsigned int * offsets, - double * out) +vectorized_transpose_and_store(const bool add_into, + const unsigned int n_entries, + const VectorizedArray *in, + const unsigned int * offsets, + double * out) { // as for the load, we split the store operations into 256 bit units to // better balance between code size, shuffle instructions, and stores @@ -996,7 +1043,7 @@ vectorized_transpose_and_store(const bool add_into, * Specialization for float and AVX512. */ template <> -class VectorizedArray +class VectorizedArray { public: /** @@ -1262,18 +1309,20 @@ private: /** * Make a few functions friends. */ - template - friend VectorizedArray - std::sqrt(const VectorizedArray &); - template - friend VectorizedArray - std::abs(const VectorizedArray &); - template - friend VectorizedArray - std::max(const VectorizedArray &, const VectorizedArray &); - template - friend VectorizedArray - std::min(const VectorizedArray &, const VectorizedArray &); + template + friend VectorizedArray + std::sqrt(const VectorizedArray &); + template + friend VectorizedArray + std::abs(const VectorizedArray &); + template + friend VectorizedArray + std::max(const VectorizedArray &, + const VectorizedArray &); + template + friend VectorizedArray + std::min(const VectorizedArray &, + const VectorizedArray &); }; @@ -1283,10 +1332,10 @@ private: */ template <> inline DEAL_II_ALWAYS_INLINE void -vectorized_load_and_transpose(const unsigned int n_entries, - const float * in, - const unsigned int * offsets, - VectorizedArray *out) +vectorized_load_and_transpose(const unsigned int n_entries, + const float * in, + const unsigned int * offsets, + VectorizedArray *out) { // Similar to the double case, we perform the work on smaller entities. In // this case, we start from 128 bit arrays and insert them into a full 512 @@ -1346,11 +1395,11 @@ vectorized_load_and_transpose(const unsigned int n_entries, */ template <> inline DEAL_II_ALWAYS_INLINE void -vectorized_transpose_and_store(const bool add_into, - const unsigned int n_entries, - const VectorizedArray *in, - const unsigned int * offsets, - float * out) +vectorized_transpose_and_store(const bool add_into, + const unsigned int n_entries, + const VectorizedArray *in, + const unsigned int * offsets, + float * out) { const unsigned int n_chunks = n_entries / 4; for (unsigned int i = 0; i < n_chunks; ++i) @@ -1453,15 +1502,15 @@ vectorized_transpose_and_store(const bool add_into, out[offsets[v] + i] = in[i][v]; } +#endif - -#elif DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 2 && defined(__AVX__) +#if DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 2 && defined(__AVX__) /** * Specialization of VectorizedArray class for double and AVX. */ template <> -class VectorizedArray +class VectorizedArray { public: /** @@ -1718,18 +1767,20 @@ private: /** * Make a few functions friends. */ - template - friend VectorizedArray - std::sqrt(const VectorizedArray &); - template - friend VectorizedArray - std::abs(const VectorizedArray &); - template - friend VectorizedArray - std::max(const VectorizedArray &, const VectorizedArray &); - template - friend VectorizedArray - std::min(const VectorizedArray &, const VectorizedArray &); + template + friend VectorizedArray + std::sqrt(const VectorizedArray &); + template + friend VectorizedArray + std::abs(const VectorizedArray &); + template + friend VectorizedArray + std::max(const VectorizedArray &, + const VectorizedArray &); + template + friend VectorizedArray + std::min(const VectorizedArray &, + const VectorizedArray &); }; @@ -1739,10 +1790,10 @@ private: */ template <> inline DEAL_II_ALWAYS_INLINE void -vectorized_load_and_transpose(const unsigned int n_entries, - const double * in, - const unsigned int * offsets, - VectorizedArray *out) +vectorized_load_and_transpose(const unsigned int n_entries, + const double * in, + const unsigned int * offsets, + VectorizedArray *out) { const unsigned int n_chunks = n_entries / 4; const double * in0 = in + offsets[0]; @@ -1778,11 +1829,11 @@ vectorized_load_and_transpose(const unsigned int n_entries, */ template <> inline DEAL_II_ALWAYS_INLINE void -vectorized_transpose_and_store(const bool add_into, - const unsigned int n_entries, - const VectorizedArray *in, - const unsigned int * offsets, - double * out) +vectorized_transpose_and_store(const bool add_into, + const unsigned int n_entries, + const VectorizedArray *in, + const unsigned int * offsets, + double * out) { const unsigned int n_chunks = n_entries / 4; double * out0 = out + offsets[0]; @@ -1844,7 +1895,7 @@ vectorized_transpose_and_store(const bool add_into, * Specialization for float and AVX. */ template <> -class VectorizedArray +class VectorizedArray { public: /** @@ -2101,18 +2152,20 @@ private: /** * Make a few functions friends. */ - template - friend VectorizedArray - std::sqrt(const VectorizedArray &); - template - friend VectorizedArray - std::abs(const VectorizedArray &); - template - friend VectorizedArray - std::max(const VectorizedArray &, const VectorizedArray &); - template - friend VectorizedArray - std::min(const VectorizedArray &, const VectorizedArray &); + template + friend VectorizedArray + std::sqrt(const VectorizedArray &); + template + friend VectorizedArray + std::abs(const VectorizedArray &); + template + friend VectorizedArray + std::max(const VectorizedArray &, + const VectorizedArray &); + template + friend VectorizedArray + std::min(const VectorizedArray &, + const VectorizedArray &); }; @@ -2122,10 +2175,10 @@ private: */ template <> inline DEAL_II_ALWAYS_INLINE void -vectorized_load_and_transpose(const unsigned int n_entries, - const float * in, - const unsigned int * offsets, - VectorizedArray *out) +vectorized_load_and_transpose(const unsigned int n_entries, + const float * in, + const unsigned int * offsets, + VectorizedArray *out) { const unsigned int n_chunks = n_entries / 4; for (unsigned int i = 0; i < n_chunks; ++i) @@ -2164,11 +2217,11 @@ vectorized_load_and_transpose(const unsigned int n_entries, */ template <> inline DEAL_II_ALWAYS_INLINE void -vectorized_transpose_and_store(const bool add_into, - const unsigned int n_entries, - const VectorizedArray *in, - const unsigned int * offsets, - float * out) +vectorized_transpose_and_store(const bool add_into, + const unsigned int n_entries, + const VectorizedArray *in, + const unsigned int * offsets, + float * out) { const unsigned int n_chunks = n_entries / 4; for (unsigned int i = 0; i < n_chunks; ++i) @@ -2240,15 +2293,15 @@ vectorized_transpose_and_store(const bool add_into, out[offsets[v] + i] = in[i][v]; } +#endif - -#elif DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 1 && defined(__SSE2__) +#if DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 1 && defined(__SSE2__) /** * Specialization for double and SSE2. */ template <> -class VectorizedArray +class VectorizedArray { public: /** @@ -2491,18 +2544,20 @@ private: /** * Make a few functions friends. */ - template - friend VectorizedArray - std::sqrt(const VectorizedArray &); - template - friend VectorizedArray - std::abs(const VectorizedArray &); - template - friend VectorizedArray - std::max(const VectorizedArray &, const VectorizedArray &); - template - friend VectorizedArray - std::min(const VectorizedArray &, const VectorizedArray &); + template + friend VectorizedArray + std::sqrt(const VectorizedArray &); + template + friend VectorizedArray + std::abs(const VectorizedArray &); + template + friend VectorizedArray + std::max(const VectorizedArray &, + const VectorizedArray &); + template + friend VectorizedArray + std::min(const VectorizedArray &, + const VectorizedArray &); }; @@ -2512,10 +2567,10 @@ private: */ template <> inline DEAL_II_ALWAYS_INLINE void -vectorized_load_and_transpose(const unsigned int n_entries, - const double * in, - const unsigned int * offsets, - VectorizedArray *out) +vectorized_load_and_transpose(const unsigned int n_entries, + const double * in, + const unsigned int * offsets, + VectorizedArray *out) { const unsigned int n_chunks = n_entries / 2; for (unsigned int i = 0; i < n_chunks; ++i) @@ -2539,11 +2594,11 @@ vectorized_load_and_transpose(const unsigned int n_entries, */ template <> inline DEAL_II_ALWAYS_INLINE void -vectorized_transpose_and_store(const bool add_into, - const unsigned int n_entries, - const VectorizedArray *in, - const unsigned int * offsets, - double * out) +vectorized_transpose_and_store(const bool add_into, + const unsigned int n_entries, + const VectorizedArray *in, + const unsigned int * offsets, + double * out) { const unsigned int n_chunks = n_entries / 2; if (add_into) @@ -2590,7 +2645,7 @@ vectorized_transpose_and_store(const bool add_into, * Specialization for float and SSE2. */ template <> -class VectorizedArray +class VectorizedArray { public: /** @@ -2833,18 +2888,20 @@ private: /** * Make a few functions friends. */ - template - friend VectorizedArray - std::sqrt(const VectorizedArray &); - template - friend VectorizedArray - std::abs(const VectorizedArray &); - template - friend VectorizedArray - std::max(const VectorizedArray &, const VectorizedArray &); - template - friend VectorizedArray - std::min(const VectorizedArray &, const VectorizedArray &); + template + friend VectorizedArray + std::sqrt(const VectorizedArray &); + template + friend VectorizedArray + std::abs(const VectorizedArray &); + template + friend VectorizedArray + std::max(const VectorizedArray &, + const VectorizedArray &); + template + friend VectorizedArray + std::min(const VectorizedArray &, + const VectorizedArray &); }; @@ -2854,10 +2911,10 @@ private: */ template <> inline DEAL_II_ALWAYS_INLINE void -vectorized_load_and_transpose(const unsigned int n_entries, - const float * in, - const unsigned int * offsets, - VectorizedArray *out) +vectorized_load_and_transpose(const unsigned int n_entries, + const float * in, + const unsigned int * offsets, + VectorizedArray *out) { const unsigned int n_chunks = n_entries / 4; for (unsigned int i = 0; i < n_chunks; ++i) @@ -2889,11 +2946,11 @@ vectorized_load_and_transpose(const unsigned int n_entries, */ template <> inline DEAL_II_ALWAYS_INLINE void -vectorized_transpose_and_store(const bool add_into, - const unsigned int n_entries, - const VectorizedArray *in, - const unsigned int * offsets, - float * out) +vectorized_transpose_and_store(const bool add_into, + const unsigned int n_entries, + const VectorizedArray *in, + const unsigned int * offsets, + float * out) { const unsigned int n_chunks = n_entries / 4; for (unsigned int i = 0; i < n_chunks; ++i) @@ -2949,12 +3006,11 @@ vectorized_transpose_and_store(const bool add_into, #endif // if DEAL_II_COMPILER_VECTORIZATION_LEVEL > 0 && defined(__SSE2__) - #if DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 1 && defined(__ALTIVEC__) && \ defined(__VSX__) template <> -class VectorizedArray +class VectorizedArray { public: /** @@ -3150,24 +3206,26 @@ private: /** * Make a few functions friends. */ - template - friend VectorizedArray - std::sqrt(const VectorizedArray &); - template - friend VectorizedArray - std::abs(const VectorizedArray &); - template - friend VectorizedArray - std::max(const VectorizedArray &, const VectorizedArray &); - template - friend VectorizedArray - std::min(const VectorizedArray &, const VectorizedArray &); + template + friend VectorizedArray + std::sqrt(const VectorizedArray &); + template + friend VectorizedArray + std::abs(const VectorizedArray &); + template + friend VectorizedArray + std::max(const VectorizedArray &, + const VectorizedArray &); + template + friend VectorizedArray + std::min(const VectorizedArray &, + const VectorizedArray &); }; template <> -class VectorizedArray +class VectorizedArray { public: /** @@ -3363,18 +3421,20 @@ private: /** * Make a few functions friends. */ - template - friend VectorizedArray - std::sqrt(const VectorizedArray &); - template - friend VectorizedArray - std::abs(const VectorizedArray &); - template - friend VectorizedArray - std::max(const VectorizedArray &, const VectorizedArray &); - template - friend VectorizedArray - std::min(const VectorizedArray &, const VectorizedArray &); + template + friend VectorizedArray + std::sqrt(const VectorizedArray &); + template + friend VectorizedArray + std::abs(const VectorizedArray &); + template + friend VectorizedArray + std::max(const VectorizedArray &, + const VectorizedArray &); + template + friend VectorizedArray + std::min(const VectorizedArray &, + const VectorizedArray &); }; #endif // if DEAL_II_VECTORIZATION_LEVEL >=1 && defined(__ALTIVEC__) && @@ -3387,12 +3447,13 @@ private: * * @relatesalso VectorizedArray */ -template +template inline DEAL_II_ALWAYS_INLINE bool -operator==(const VectorizedArray &lhs, - const VectorizedArray &rhs) +operator==(const VectorizedArray &lhs, + const VectorizedArray &rhs) { - for (unsigned int i = 0; i < VectorizedArray::n_array_elements; ++i) + for (unsigned int i = 0; i < VectorizedArray::n_array_elements; + ++i) if (lhs[i] != rhs[i]) return false; @@ -3405,11 +3466,12 @@ operator==(const VectorizedArray &lhs, * * @relatesalso VectorizedArray */ -template -inline DEAL_II_ALWAYS_INLINE VectorizedArray - operator+(const VectorizedArray &u, const VectorizedArray &v) +template +inline DEAL_II_ALWAYS_INLINE VectorizedArray + operator+(const VectorizedArray &u, + const VectorizedArray &v) { - VectorizedArray tmp = u; + VectorizedArray tmp = u; return tmp += v; } @@ -3418,11 +3480,12 @@ inline DEAL_II_ALWAYS_INLINE VectorizedArray * * @relatesalso VectorizedArray */ -template -inline DEAL_II_ALWAYS_INLINE VectorizedArray - operator-(const VectorizedArray &u, const VectorizedArray &v) +template +inline DEAL_II_ALWAYS_INLINE VectorizedArray + operator-(const VectorizedArray &u, + const VectorizedArray &v) { - VectorizedArray tmp = u; + VectorizedArray tmp = u; return tmp -= v; } @@ -3431,11 +3494,12 @@ inline DEAL_II_ALWAYS_INLINE VectorizedArray * * @relatesalso VectorizedArray */ -template -inline DEAL_II_ALWAYS_INLINE VectorizedArray - operator*(const VectorizedArray &u, const VectorizedArray &v) +template +inline DEAL_II_ALWAYS_INLINE VectorizedArray + operator*(const VectorizedArray &u, + const VectorizedArray &v) { - VectorizedArray tmp = u; + VectorizedArray tmp = u; return tmp *= v; } @@ -3444,11 +3508,12 @@ inline DEAL_II_ALWAYS_INLINE VectorizedArray * * @relatesalso VectorizedArray */ -template -inline DEAL_II_ALWAYS_INLINE VectorizedArray - operator/(const VectorizedArray &u, const VectorizedArray &v) +template +inline DEAL_II_ALWAYS_INLINE VectorizedArray + operator/(const VectorizedArray &u, + const VectorizedArray &v) { - VectorizedArray tmp = u; + VectorizedArray tmp = u; return tmp /= v; } @@ -3458,11 +3523,11 @@ inline DEAL_II_ALWAYS_INLINE VectorizedArray * * @relatesalso VectorizedArray */ -template -inline DEAL_II_ALWAYS_INLINE VectorizedArray - operator+(const Number &u, const VectorizedArray &v) +template +inline DEAL_II_ALWAYS_INLINE VectorizedArray + operator+(const Number &u, const VectorizedArray &v) { - VectorizedArray tmp; + VectorizedArray tmp; tmp = u; return tmp += v; } @@ -3475,10 +3540,11 @@ inline DEAL_II_ALWAYS_INLINE VectorizedArray * * @relatesalso VectorizedArray */ -inline DEAL_II_ALWAYS_INLINE VectorizedArray - operator+(const double u, const VectorizedArray &v) +template +inline DEAL_II_ALWAYS_INLINE VectorizedArray + operator+(const double u, const VectorizedArray &v) { - VectorizedArray tmp; + VectorizedArray tmp; tmp = u; return tmp += v; } @@ -3489,9 +3555,9 @@ inline DEAL_II_ALWAYS_INLINE VectorizedArray * * @relatesalso VectorizedArray */ -template -inline DEAL_II_ALWAYS_INLINE VectorizedArray - operator+(const VectorizedArray &v, const Number &u) +template +inline DEAL_II_ALWAYS_INLINE VectorizedArray + operator+(const VectorizedArray &v, const Number &u) { return u + v; } @@ -3504,8 +3570,9 @@ inline DEAL_II_ALWAYS_INLINE VectorizedArray * * @relatesalso VectorizedArray */ -inline DEAL_II_ALWAYS_INLINE VectorizedArray - operator+(const VectorizedArray &v, const double u) +template +inline DEAL_II_ALWAYS_INLINE VectorizedArray + operator+(const VectorizedArray &v, const double u) { return u + v; } @@ -3516,11 +3583,11 @@ inline DEAL_II_ALWAYS_INLINE VectorizedArray * * @relatesalso VectorizedArray */ -template -inline DEAL_II_ALWAYS_INLINE VectorizedArray - operator-(const Number &u, const VectorizedArray &v) +template +inline DEAL_II_ALWAYS_INLINE VectorizedArray + operator-(const Number &u, const VectorizedArray &v) { - VectorizedArray tmp; + VectorizedArray tmp; tmp = u; return tmp -= v; } @@ -3533,10 +3600,11 @@ inline DEAL_II_ALWAYS_INLINE VectorizedArray * * @relatesalso VectorizedArray */ -inline DEAL_II_ALWAYS_INLINE VectorizedArray - operator-(const double u, const VectorizedArray &v) +template +inline DEAL_II_ALWAYS_INLINE VectorizedArray + operator-(const double u, const VectorizedArray &v) { - VectorizedArray tmp; + VectorizedArray tmp; tmp = float(u); return tmp -= v; } @@ -3547,11 +3615,11 @@ inline DEAL_II_ALWAYS_INLINE VectorizedArray * * @relatesalso VectorizedArray */ -template -inline DEAL_II_ALWAYS_INLINE VectorizedArray - operator-(const VectorizedArray &v, const Number &u) +template +inline DEAL_II_ALWAYS_INLINE VectorizedArray + operator-(const VectorizedArray &v, const Number &u) { - VectorizedArray tmp; + VectorizedArray tmp; tmp = u; return v - tmp; } @@ -3564,10 +3632,11 @@ inline DEAL_II_ALWAYS_INLINE VectorizedArray * * @relatesalso VectorizedArray */ -inline DEAL_II_ALWAYS_INLINE VectorizedArray - operator-(const VectorizedArray &v, const double u) +template +inline DEAL_II_ALWAYS_INLINE VectorizedArray + operator-(const VectorizedArray &v, const double u) { - VectorizedArray tmp; + VectorizedArray tmp; tmp = float(u); return v - tmp; } @@ -3578,11 +3647,11 @@ inline DEAL_II_ALWAYS_INLINE VectorizedArray * * @relatesalso VectorizedArray */ -template -inline DEAL_II_ALWAYS_INLINE VectorizedArray - operator*(const Number &u, const VectorizedArray &v) +template +inline DEAL_II_ALWAYS_INLINE VectorizedArray + operator*(const Number &u, const VectorizedArray &v) { - VectorizedArray tmp; + VectorizedArray tmp; tmp = u; return tmp *= v; } @@ -3595,10 +3664,11 @@ inline DEAL_II_ALWAYS_INLINE VectorizedArray * * @relatesalso VectorizedArray */ -inline DEAL_II_ALWAYS_INLINE VectorizedArray - operator*(const double u, const VectorizedArray &v) +template +inline DEAL_II_ALWAYS_INLINE VectorizedArray + operator*(const double u, const VectorizedArray &v) { - VectorizedArray tmp; + VectorizedArray tmp; tmp = float(u); return tmp *= v; } @@ -3609,9 +3679,9 @@ inline DEAL_II_ALWAYS_INLINE VectorizedArray * * @relatesalso VectorizedArray */ -template -inline DEAL_II_ALWAYS_INLINE VectorizedArray - operator*(const VectorizedArray &v, const Number &u) +template +inline DEAL_II_ALWAYS_INLINE VectorizedArray + operator*(const VectorizedArray &v, const Number &u) { return u * v; } @@ -3624,8 +3694,9 @@ inline DEAL_II_ALWAYS_INLINE VectorizedArray * * @relatesalso VectorizedArray */ -inline DEAL_II_ALWAYS_INLINE VectorizedArray - operator*(const VectorizedArray &v, const double u) +template +inline DEAL_II_ALWAYS_INLINE VectorizedArray + operator*(const VectorizedArray &v, const double u) { return u * v; } @@ -3636,11 +3707,11 @@ inline DEAL_II_ALWAYS_INLINE VectorizedArray * * @relatesalso VectorizedArray */ -template -inline DEAL_II_ALWAYS_INLINE VectorizedArray - operator/(const Number &u, const VectorizedArray &v) +template +inline DEAL_II_ALWAYS_INLINE VectorizedArray + operator/(const Number &u, const VectorizedArray &v) { - VectorizedArray tmp; + VectorizedArray tmp; tmp = u; return tmp /= v; } @@ -3653,10 +3724,11 @@ inline DEAL_II_ALWAYS_INLINE VectorizedArray * * @relatesalso VectorizedArray */ -inline DEAL_II_ALWAYS_INLINE VectorizedArray - operator/(const double u, const VectorizedArray &v) +template +inline DEAL_II_ALWAYS_INLINE VectorizedArray + operator/(const double u, const VectorizedArray &v) { - VectorizedArray tmp; + VectorizedArray tmp; tmp = float(u); return tmp /= v; } @@ -3667,11 +3739,11 @@ inline DEAL_II_ALWAYS_INLINE VectorizedArray * * @relatesalso VectorizedArray */ -template -inline DEAL_II_ALWAYS_INLINE VectorizedArray - operator/(const VectorizedArray &v, const Number &u) +template +inline DEAL_II_ALWAYS_INLINE VectorizedArray + operator/(const VectorizedArray &v, const Number &u) { - VectorizedArray tmp; + VectorizedArray tmp; tmp = u; return v / tmp; } @@ -3684,10 +3756,11 @@ inline DEAL_II_ALWAYS_INLINE VectorizedArray * * @relatesalso VectorizedArray */ -inline DEAL_II_ALWAYS_INLINE VectorizedArray - operator/(const VectorizedArray &v, const double u) +template +inline DEAL_II_ALWAYS_INLINE VectorizedArray + operator/(const VectorizedArray &v, const double u) { - VectorizedArray tmp; + VectorizedArray tmp; tmp = float(u); return v / tmp; } @@ -3697,9 +3770,9 @@ inline DEAL_II_ALWAYS_INLINE VectorizedArray * * @relatesalso VectorizedArray */ -template -inline DEAL_II_ALWAYS_INLINE VectorizedArray - operator+(const VectorizedArray &u) +template +inline DEAL_II_ALWAYS_INLINE VectorizedArray + operator+(const VectorizedArray &u) { return u; } @@ -3709,13 +3782,13 @@ inline DEAL_II_ALWAYS_INLINE VectorizedArray * * @relatesalso VectorizedArray */ -template -inline DEAL_II_ALWAYS_INLINE VectorizedArray - operator-(const VectorizedArray &u) +template +inline DEAL_II_ALWAYS_INLINE VectorizedArray + operator-(const VectorizedArray &u) { // to get a negative sign, subtract the input from zero (could also // multiply by -1, but this one is slightly simpler) - return VectorizedArray() - u; + return VectorizedArray() - u; } @@ -3737,21 +3810,21 @@ namespace std * * @relatesalso VectorizedArray */ - template - inline ::dealii::VectorizedArray - sin(const ::dealii::VectorizedArray &x) + template + inline ::dealii::VectorizedArray + sin(const ::dealii::VectorizedArray &x) { // put values in an array and later read in that array with an unaligned // 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 with SSE2 (see also deal.II developers list // from April 2014, topic "matrix_free/step-48 Test"). - Number values[::dealii::VectorizedArray::n_array_elements]; + Number values[::dealii::VectorizedArray::n_array_elements]; for (unsigned int i = 0; - i < dealii::VectorizedArray::n_array_elements; + i < dealii::VectorizedArray::n_array_elements; ++i) values[i] = std::sin(x[i]); - ::dealii::VectorizedArray out; + ::dealii::VectorizedArray out; out.load(&values[0]); return out; } @@ -3765,16 +3838,16 @@ namespace std * * @relatesalso VectorizedArray */ - template - inline ::dealii::VectorizedArray - cos(const ::dealii::VectorizedArray &x) + template + inline ::dealii::VectorizedArray + cos(const ::dealii::VectorizedArray &x) { - Number values[::dealii::VectorizedArray::n_array_elements]; + Number values[::dealii::VectorizedArray::n_array_elements]; for (unsigned int i = 0; - i < dealii::VectorizedArray::n_array_elements; + i < dealii::VectorizedArray::n_array_elements; ++i) values[i] = std::cos(x[i]); - ::dealii::VectorizedArray out; + ::dealii::VectorizedArray out; out.load(&values[0]); return out; } @@ -3788,16 +3861,16 @@ namespace std * * @relatesalso VectorizedArray */ - template - inline ::dealii::VectorizedArray - tan(const ::dealii::VectorizedArray &x) + template + inline ::dealii::VectorizedArray + tan(const ::dealii::VectorizedArray &x) { - Number values[::dealii::VectorizedArray::n_array_elements]; + Number values[::dealii::VectorizedArray::n_array_elements]; for (unsigned int i = 0; - i < dealii::VectorizedArray::n_array_elements; + i < dealii::VectorizedArray::n_array_elements; ++i) values[i] = std::tan(x[i]); - ::dealii::VectorizedArray out; + ::dealii::VectorizedArray out; out.load(&values[0]); return out; } @@ -3811,16 +3884,16 @@ namespace std * * @relatesalso VectorizedArray */ - template - inline ::dealii::VectorizedArray - exp(const ::dealii::VectorizedArray &x) + template + inline ::dealii::VectorizedArray + exp(const ::dealii::VectorizedArray &x) { - Number values[::dealii::VectorizedArray::n_array_elements]; + Number values[::dealii::VectorizedArray::n_array_elements]; for (unsigned int i = 0; - i < dealii::VectorizedArray::n_array_elements; + i < dealii::VectorizedArray::n_array_elements; ++i) values[i] = std::exp(x[i]); - ::dealii::VectorizedArray out; + ::dealii::VectorizedArray out; out.load(&values[0]); return out; } @@ -3834,16 +3907,16 @@ namespace std * * @relatesalso VectorizedArray */ - template - inline ::dealii::VectorizedArray - log(const ::dealii::VectorizedArray &x) + template + inline ::dealii::VectorizedArray + log(const ::dealii::VectorizedArray &x) { - Number values[::dealii::VectorizedArray::n_array_elements]; + Number values[::dealii::VectorizedArray::n_array_elements]; for (unsigned int i = 0; - i < dealii::VectorizedArray::n_array_elements; + i < dealii::VectorizedArray::n_array_elements; ++i) values[i] = std::log(x[i]); - ::dealii::VectorizedArray out; + ::dealii::VectorizedArray out; out.load(&values[0]); return out; } @@ -3857,9 +3930,9 @@ namespace std * * @relatesalso VectorizedArray */ - template - inline ::dealii::VectorizedArray - sqrt(const ::dealii::VectorizedArray &x) + template + inline ::dealii::VectorizedArray + sqrt(const ::dealii::VectorizedArray &x) { return x.get_sqrt(); } @@ -3873,16 +3946,16 @@ namespace std * * @relatesalso VectorizedArray */ - template - inline ::dealii::VectorizedArray - pow(const ::dealii::VectorizedArray &x, const Number p) + template + inline ::dealii::VectorizedArray + pow(const ::dealii::VectorizedArray &x, const Number p) { - Number values[::dealii::VectorizedArray::n_array_elements]; + Number values[::dealii::VectorizedArray::n_array_elements]; for (unsigned int i = 0; - i < dealii::VectorizedArray::n_array_elements; + i < dealii::VectorizedArray::n_array_elements; ++i) values[i] = std::pow(x[i], p); - ::dealii::VectorizedArray out; + ::dealii::VectorizedArray out; out.load(&values[0]); return out; } @@ -3896,9 +3969,9 @@ namespace std * * @relatesalso VectorizedArray */ - template - inline ::dealii::VectorizedArray - abs(const ::dealii::VectorizedArray &x) + template + inline ::dealii::VectorizedArray + abs(const ::dealii::VectorizedArray &x) { return x.get_abs(); } @@ -3912,10 +3985,10 @@ namespace std * * @relatesalso VectorizedArray */ - template - inline ::dealii::VectorizedArray - max(const ::dealii::VectorizedArray &x, - const ::dealii::VectorizedArray &y) + template + inline ::dealii::VectorizedArray + max(const ::dealii::VectorizedArray &x, + const ::dealii::VectorizedArray &y) { return x.get_max(y); } @@ -3929,10 +4002,10 @@ namespace std * * @relatesalso VectorizedArray */ - template - inline ::dealii::VectorizedArray - min(const ::dealii::VectorizedArray &x, - const ::dealii::VectorizedArray &y) + template + inline ::dealii::VectorizedArray + min(const ::dealii::VectorizedArray &x, + const ::dealii::VectorizedArray &y) { return x.get_min(y); } diff --git a/include/deal.II/lac/tensor_product_matrix.h b/include/deal.II/lac/tensor_product_matrix.h index 3e5a671515..b6b1f9e108 100644 --- a/include/deal.II/lac/tensor_product_matrix.h +++ b/include/deal.II/lac/tensor_product_matrix.h @@ -32,8 +32,6 @@ template class Vector; template class FullMatrix; -template -class VectorizedArray; /** * This is an abstract base class used for a special matrix class, namely the diff --git a/source/base/vectorization.cc b/source/base/vectorization.cc index 5d544e15d8..f05ecb7b3a 100644 --- a/source/base/vectorization.cc +++ b/source/base/vectorization.cc @@ -17,11 +17,6 @@ DEAL_II_NAMESPACE_OPEN -#if DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 1 && !defined(DEAL_II_MSVC) -const unsigned int VectorizedArray::n_array_elements; -const unsigned int VectorizedArray::n_array_elements; -#endif - // VectorizedArray must be a POD (plain old data) type to make sure it // can use maximum level of compiler optimization. // A type is POD if it has standard layout (similar to a C struct) @@ -35,4 +30,33 @@ static_assert(std::is_standard_layout>::value && std::is_trivial>::value, "VectorizedArray must be a POD type"); +#if DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 1 && !defined(DEAL_II_MSVC) +# if DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 3 && defined(__AVX512F__) +const unsigned int VectorizedArray::n_array_elements; +# endif + +# if DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 2 && defined(__AVX__) +const unsigned int VectorizedArray::n_array_elements; +# endif + +# if (DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 1 && defined(__SSE2__)) || \ + (DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 1 && defined(__ALTIVEC__)) +const unsigned int VectorizedArray::n_array_elements; +# endif + + +# if DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 3 && defined(__AVX512F__) +const unsigned int VectorizedArray::n_array_elements; +# endif + +# if DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 2 && defined(__AVX__) +const unsigned int VectorizedArray::n_array_elements; +# endif + +# if (DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 1 && defined(__SSE2__)) || \ + (DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 1 && defined(__ALTIVEC__)) +const unsigned int VectorizedArray::n_array_elements; +# endif +#endif + DEAL_II_NAMESPACE_CLOSE diff --git a/tests/base/vectorization_10.cc b/tests/base/vectorization_10.cc new file mode 100644 index 0000000000..6823059429 --- /dev/null +++ b/tests/base/vectorization_10.cc @@ -0,0 +1,182 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2015 - 2018 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + +// test transpose operations of vectorized array using the array+offset method +// for the set of all supported vectorization widths (otherwise the same as +// vectorization_05) + +#include + +#include + +#include "../tests.h" + + +template +void +do_test() +{ + // since the number of array elements is system dependent, it is not a good + // idea to print them to an output file. Instead, check the values manually + const unsigned int n_vectors = + VectorizedArray::n_array_elements; + VectorizedArray arr[n_numbers]; + Number other[n_vectors * n_numbers]; + unsigned int offsets[n_vectors]; + for (unsigned int v = 0; v < n_vectors; ++v) + offsets[v] = v * n_numbers; + + for (unsigned int i = 0; i < n_vectors; ++i) + for (unsigned int j = 0; j < n_numbers; ++j) + other[i * n_numbers + j] = i * n_numbers + j; + + vectorized_load_and_transpose(n_numbers, other, offsets, arr); + unsigned int n_errors = 0; + for (unsigned int j = 0; j < n_numbers; ++j) + for (unsigned int i = 0; i < n_vectors; ++i) + if (arr[j][i] != i * n_numbers + j) + ++n_errors; + if (n_errors > 0) + { + deallog << "load_and_transpose at n=" << n_numbers + << " width=" << width << ": #errors: " << n_errors << std::endl; + + for (unsigned int i = 0; i < n_numbers; ++i) + { + for (unsigned int j = 0; j < n_vectors; ++j) + deallog << arr[i][j] << " "; + deallog << std::endl; + } + } + + vectorized_transpose_and_store(true, n_numbers, arr, offsets, other); + n_errors = 0; + for (unsigned int i = 0; i < n_vectors; ++i) + for (unsigned int j = 0; j < n_numbers; ++j) + if (other[i * n_numbers + j] != 2. * (i * n_numbers + j)) + ++n_errors; + if (n_errors > 0) + { + deallog << "transpose_and_store ( add) at n=" << n_numbers + << " width=" << width << ": #errors: " << n_errors << std::endl; + + for (unsigned int i = 0; i < n_vectors; ++i) + { + for (unsigned int j = 0; j < n_numbers; ++j) + deallog << other[i * n_numbers + j] << " "; + deallog << std::endl; + } + } + + vectorized_transpose_and_store(false, n_numbers, arr, offsets, other); + n_errors = 0; + for (unsigned int i = 0; i < n_vectors; ++i) + for (unsigned int j = 0; j < n_numbers; ++j) + if (other[i * n_numbers + j] != (i * n_numbers + j)) + ++n_errors; + if (n_errors > 0) + { + deallog << "transpose_and_store (noadd) at n=" << n_numbers + << " width=" << width << ": #errors: " << n_errors << std::endl; + + for (unsigned int i = 0; i < n_vectors; ++i) + { + for (unsigned int j = 0; j < n_numbers; ++j) + deallog << other[i * n_numbers + j] << " "; + deallog << std::endl; + } + } +} + + +template +struct Tester +{ + static void + test() + { + do_test(); + } +}; + +template +struct Tester +{ + static void + test() + { +#if DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 3 && defined(__AVX512F__) + do_test(); +#endif + +#if DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 2 && defined(__AVX__) + do_test(); +#endif + +#if DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 1 && defined(__SSE2__) + do_test(); +#endif + + do_test(); + } +}; + +template +struct Tester +{ + static void + test() + { +#if DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 3 && defined(__AVX512F__) + do_test(); +#endif + +#if DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 2 && defined(__AVX__) + do_test(); +#endif + +#if DEAL_II_COMPILER_VECTORIZATION_LEVEL >= 1 && defined(__SSE2__) + do_test(); +#endif + + do_test(); + } +}; + + +int +main() +{ + initlog(); + + deallog.push("double"); + Tester::test(); + Tester::test(); + Tester::test(); + deallog.pop(); + deallog.push("float"); + Tester::test(); + Tester::test(); + Tester::test(); + deallog.pop(); + + // test long double: in that case, the default + // path of VectorizedArray is taken no matter + // what was done for double or float + deallog.push("long double"); + Tester::test(); + deallog.pop(); +} diff --git a/tests/base/vectorization_10.output b/tests/base/vectorization_10.output new file mode 100644 index 0000000000..8b13789179 --- /dev/null +++ b/tests/base/vectorization_10.output @@ -0,0 +1 @@ +