MCKL
Monte Carlo Kernel Library
threefry_avx2_2x32.hpp
Go to the documentation of this file.
1 //============================================================================
2 // MCKL/include/mckl/random/internal/threefry_avx2_2x32.hpp
3 //----------------------------------------------------------------------------
4 // MCKL: Monte Carlo Kernel Library
5 //----------------------------------------------------------------------------
6 // Copyright (c) 2013-2018, Yan Zhou
7 // All rights reserved.
8 //
9 // Redistribution and use in source and binary forms, with or without
10 // modification, are permitted provided that the following conditions are met:
11 //
12 // Redistributions of source code must retain the above copyright notice,
13 // this list of conditions and the following disclaimer.
14 //
15 // Redistributions in binary form must reproduce the above copyright notice,
16 // this list of conditions and the following disclaimer in the documentation
17 // and/or other materials provided with the distribution.
18 //
19 // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS AS IS
20 // AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
21 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
22 // ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
23 // LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
24 // CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
25 // SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
26 // INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
27 // CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
28 // ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
29 // POSSIBILITY OF SUCH DAMAGE.
30 //============================================================================
31 
32 #ifndef MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_HPP
33 #define MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_HPP
34 
35 #include <mckl/random/internal/common.hpp>
36 #include <mckl/random/internal/threefry_common.hpp>
37 #include <mckl/random/internal/threefry_constants.hpp>
38 #include <mckl/random/internal/threefry_generic_2x32.hpp>
39 
40 MCKL_PUSH_GCC_WARNING("-Wignored-attributes")
41 
42 #define MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_KBOX(N) \
43  ymmt0 = _mm256_set1_epi32( \
44  static_cast<int>(ThreefryKBox<T, K, N>::template key<0>(par))); \
45  ymmt1 = _mm256_set1_epi32( \
46  static_cast<int>(ThreefryKBox<T, K, N>::template key<1>(par))); \
47  ymms0 = _mm256_add_epi32(ymms0, ymmt0); \
48  ymms1 = _mm256_add_epi32(ymms1, ymmt1); \
49  ymms2 = _mm256_add_epi32(ymms2, ymmt0); \
50  ymms3 = _mm256_add_epi32(ymms3, ymmt1); \
51  ymms4 = _mm256_add_epi32(ymms4, ymmt0); \
52  ymms5 = _mm256_add_epi32(ymms5, ymmt1); \
53  ymms6 = _mm256_add_epi32(ymms6, ymmt0); \
54  ymms7 = _mm256_add_epi32(ymms7, ymmt1);
55 
56 #define MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_RBOX(N) \
57  { \
58  constexpr int L0 = Constants::rotate::value[0][(N - 1) % 8]; \
59  constexpr int R0 = 32 - L0; \
60  ymms0 = _mm256_add_epi32(ymms0, ymms1); \
61  ymms2 = _mm256_add_epi32(ymms2, ymms3); \
62  ymms4 = _mm256_add_epi32(ymms4, ymms5); \
63  ymms6 = _mm256_add_epi32(ymms6, ymms7); \
64  ymmt1 = _mm256_slli_epi32(ymms1, L0); \
65  ymmt3 = _mm256_slli_epi32(ymms3, L0); \
66  ymmt5 = _mm256_slli_epi32(ymms5, L0); \
67  ymmt7 = _mm256_slli_epi32(ymms7, L0); \
68  ymms1 = _mm256_srli_epi32(ymms1, R0); \
69  ymms3 = _mm256_srli_epi32(ymms3, R0); \
70  ymms5 = _mm256_srli_epi32(ymms5, R0); \
71  ymms7 = _mm256_srli_epi32(ymms7, R0); \
72  ymmt1 = _mm256_or_si256(ymms1, ymmt1); \
73  ymmt3 = _mm256_or_si256(ymms3, ymmt3); \
74  ymmt5 = _mm256_or_si256(ymms5, ymmt5); \
75  ymmt7 = _mm256_or_si256(ymms7, ymmt7); \
76  ymms1 = _mm256_xor_si256(ymms0, ymmt1); \
77  ymms3 = _mm256_xor_si256(ymms2, ymmt3); \
78  ymms5 = _mm256_xor_si256(ymms4, ymmt5); \
79  ymms7 = _mm256_xor_si256(ymms6, ymmt7); \
80  }
81 
82 #define MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_CYCLE_4(N) \
83  MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_RBOX(N * 8 + 1); \
84  MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_RBOX(N * 8 + 2); \
85  MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_RBOX(N * 8 + 3); \
86  MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_RBOX(N * 8 + 4); \
87  MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_KBOX(N * 8 + 4);
88 
89 #define MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_CYCLE_8(N) \
90  MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_RBOX(N * 8 + 1); \
91  MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_RBOX(N * 8 + 2); \
92  MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_RBOX(N * 8 + 3); \
93  MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_RBOX(N * 8 + 4); \
94  MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_KBOX(N * 8 + 4); \
95  MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_RBOX(N * 8 + 5); \
96  MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_RBOX(N * 8 + 6); \
97  MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_RBOX(N * 8 + 7); \
98  MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_RBOX(N * 8 + 8); \
99  MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_KBOX(N * 8 + 8);
100 
101 namespace mckl {
102 
103 namespace internal {
104 
105 template <typename T>
106 class Threefry2x32GeneratorAVX2Impl
107 {
108  static_assert(std::numeric_limits<T>::digits == 32,
109  "**Threefry2x32GeneratorAVX2Impl** used with T other than a "
110  "32-bit unsigned integers");
111 
112  static constexpr std::size_t K = 2;
113  static constexpr std::size_t Rounds = 20;
114  using Constants = ThreefryConstants<T, 2>;
115 
116  public:
117  static void eval(
118  const void *plain, void *cipher, const std::array<T, K + 4> &par)
119  {
120  Threefry2x32GeneratorGenericImpl<T>::eval(plain, cipher, par);
121  }
122 
123  template <typename ResultType>
124  static void eval(std::array<std::uint64_t, 1> &ctr, ResultType *r,
125  const std::array<T, K + 4> &par)
126  {
127  Threefry2x32GeneratorGenericImpl<T>::eval(ctr, r, par);
128  }
129 
130  template <typename ResultType>
131  static void eval(std::array<std::uint64_t, 1> &ctr, std::size_t n,
132  ResultType *r, const std::array<T, K + 4> &par)
133  {
134  constexpr std::size_t R = sizeof(T) * K / sizeof(ResultType);
135 
136  const std::size_t n0 =
137  static_cast<std::size_t>(std::min(static_cast<std::uint64_t>(n),
138  std::numeric_limits<std::uint64_t>::max() - ctr.front()));
139 
140  eval_kernel(ctr, n0, r, par);
141  n -= n0;
142  r += n0 * R;
143 
144  if (n != 0) {
145  eval(ctr, r, par);
146  n -= 1;
147  r += R;
148  }
149 
150  eval_kernel(ctr, n, r, par);
151  }
152 
153  private:
154  template <typename ResultType>
155  static void eval_kernel(std::array<std::uint64_t, 1> &ctr, std::size_t n,
156  ResultType *r, const std::array<T, K + 4> &par)
157  {
158  constexpr std::size_t S = 8;
159  constexpr std::size_t N = sizeof(__m256i) * S / (sizeof(T) * K);
160 
161  __m256i ymmc =
162  _mm256_set1_epi64x(static_cast<MCKL_INT64>(std::get<0>(ctr)));
163  ctr.front() += n;
164 
165  __m256i *rptr = reinterpret_cast<__m256i *>(r);
166  while (n != 0) {
167  __m256i ymms0 = _mm256_add_epi64(
168  ymmc, _mm256_set_epi64x(0x04, 0x03, 0x02, 0x01));
169  __m256i ymms1 = _mm256_add_epi64(
170  ymmc, _mm256_set_epi64x(0x08, 0x07, 0x06, 0x05));
171  __m256i ymms2 = _mm256_add_epi64(
172  ymmc, _mm256_set_epi64x(0x0C, 0x0B, 0x0A, 0x09));
173  __m256i ymms3 = _mm256_add_epi64(
174  ymmc, _mm256_set_epi64x(0x10, 0x0F, 0x0E, 0x0D));
175  __m256i ymms4 = _mm256_add_epi64(
176  ymmc, _mm256_set_epi64x(0x14, 0x13, 0x12, 0x11));
177  __m256i ymms5 = _mm256_add_epi64(
178  ymmc, _mm256_set_epi64x(0x18, 0x17, 0x16, 0x15));
179  __m256i ymms6 = _mm256_add_epi64(
180  ymmc, _mm256_set_epi64x(0x1C, 0x1B, 0x1A, 0x19));
181  __m256i ymms7 = _mm256_add_epi64(
182  ymmc, _mm256_set_epi64x(0x20, 0x1F, 0x1E, 0x1D));
183  ymmc = _mm256_add_epi64(ymmc, _mm256_set1_epi64x(0x20));
184 
185  __m256i ymmt0 = _mm256_unpacklo_epi32(ymms0, ymms1);
186  __m256i ymmt1 = _mm256_unpacklo_epi32(ymms2, ymms3);
187  __m256i ymmt2 = _mm256_unpacklo_epi32(ymms4, ymms5);
188  __m256i ymmt3 = _mm256_unpacklo_epi32(ymms6, ymms7);
189  __m256i ymmt4 = _mm256_unpackhi_epi32(ymms0, ymms1);
190  __m256i ymmt5 = _mm256_unpackhi_epi32(ymms2, ymms3);
191  __m256i ymmt6 = _mm256_unpackhi_epi32(ymms4, ymms5);
192  __m256i ymmt7 = _mm256_unpackhi_epi32(ymms6, ymms7);
193  __m256i ymmt8 = _mm256_unpacklo_epi64(ymmt0, ymmt1);
194  __m256i ymmt9 = _mm256_unpacklo_epi64(ymmt2, ymmt3);
195  __m256i ymmtA = _mm256_unpacklo_epi64(ymmt4, ymmt5);
196  __m256i ymmtB = _mm256_unpacklo_epi64(ymmt6, ymmt7);
197  __m256i ymmtC = _mm256_unpackhi_epi64(ymmt0, ymmt1);
198  __m256i ymmtD = _mm256_unpackhi_epi64(ymmt2, ymmt3);
199  __m256i ymmtE = _mm256_unpackhi_epi64(ymmt4, ymmt5);
200  __m256i ymmtF = _mm256_unpackhi_epi64(ymmt6, ymmt7);
201  ymms0 = _mm256_permute2x128_si256(ymmt8, ymmt9, 0x20);
202  ymms1 = _mm256_permute2x128_si256(ymmtC, ymmtD, 0x20);
203  ymms2 = _mm256_permute2x128_si256(ymmtA, ymmtB, 0x20);
204  ymms3 = _mm256_permute2x128_si256(ymmtE, ymmtF, 0x20);
205  ymms4 = _mm256_permute2x128_si256(ymmt8, ymmt9, 0x31);
206  ymms5 = _mm256_permute2x128_si256(ymmtC, ymmtD, 0x31);
207  ymms6 = _mm256_permute2x128_si256(ymmtA, ymmtB, 0x31);
208  ymms7 = _mm256_permute2x128_si256(ymmtE, ymmtF, 0x31);
209 
210  MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_KBOX(0)
211  MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_CYCLE_8(0)
212  MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_CYCLE_8(1)
213  MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_CYCLE_4(2)
214 
215  ymmt0 = _mm256_unpacklo_epi32(ymms0, ymms1);
216  ymmt1 = _mm256_unpacklo_epi32(ymms2, ymms3);
217  ymmt2 = _mm256_unpacklo_epi32(ymms4, ymms5);
218  ymmt3 = _mm256_unpacklo_epi32(ymms6, ymms7);
219  ymmt4 = _mm256_unpackhi_epi32(ymms0, ymms1);
220  ymmt5 = _mm256_unpackhi_epi32(ymms2, ymms3);
221  ymmt6 = _mm256_unpackhi_epi32(ymms4, ymms5);
222  ymmt7 = _mm256_unpackhi_epi32(ymms6, ymms7);
223  ymmt8 = _mm256_unpacklo_epi64(ymmt0, ymmt1);
224  ymmt9 = _mm256_unpacklo_epi64(ymmt2, ymmt3);
225  ymmtA = _mm256_unpacklo_epi64(ymmt4, ymmt5);
226  ymmtB = _mm256_unpacklo_epi64(ymmt6, ymmt7);
227  ymmtC = _mm256_unpackhi_epi64(ymmt0, ymmt1);
228  ymmtD = _mm256_unpackhi_epi64(ymmt2, ymmt3);
229  ymmtE = _mm256_unpackhi_epi64(ymmt4, ymmt5);
230  ymmtF = _mm256_unpackhi_epi64(ymmt6, ymmt7);
231  ymms0 = _mm256_permute2x128_si256(ymmt8, ymmt9, 0x20);
232  ymms1 = _mm256_permute2x128_si256(ymmtC, ymmtD, 0x20);
233  ymms2 = _mm256_permute2x128_si256(ymmtA, ymmtB, 0x20);
234  ymms3 = _mm256_permute2x128_si256(ymmtE, ymmtF, 0x20);
235  ymms4 = _mm256_permute2x128_si256(ymmt8, ymmt9, 0x31);
236  ymms5 = _mm256_permute2x128_si256(ymmtC, ymmtD, 0x31);
237  ymms6 = _mm256_permute2x128_si256(ymmtA, ymmtB, 0x31);
238  ymms7 = _mm256_permute2x128_si256(ymmtE, ymmtF, 0x31);
239 
240  if (n >= N) {
241  n -= N;
242  _mm256_storeu_si256(rptr++, ymms0);
243  _mm256_storeu_si256(rptr++, ymms1);
244  _mm256_storeu_si256(rptr++, ymms2);
245  _mm256_storeu_si256(rptr++, ymms3);
246  _mm256_storeu_si256(rptr++, ymms4);
247  _mm256_storeu_si256(rptr++, ymms5);
248  _mm256_storeu_si256(rptr++, ymms6);
249  _mm256_storeu_si256(rptr++, ymms7);
250  } else {
251  std::array<__m256i, S> s;
252  std::get<0>(s) = ymms0;
253  std::get<1>(s) = ymms1;
254  std::get<2>(s) = ymms2;
255  std::get<3>(s) = ymms3;
256  std::get<4>(s) = ymms4;
257  std::get<5>(s) = ymms5;
258  std::get<6>(s) = ymms6;
259  std::get<7>(s) = ymms7;
260  std::memcpy(rptr, s.data(), n * sizeof(T) * K);
261  break;
262  }
263  }
264  }
265 }; // class Threefry2x32GeneratorAVX2Impl
266 
267 } // namespace internal
268 
269 } // namespace mckl
270 
271 MCKL_POP_GCC_WARNING
272 
273 #endif // MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_HPP
threefry_generic_2x32.hpp
mckl::internal::Threefry2x32GeneratorAVX2Impl::eval
static void eval(std::array< std::uint64_t, 1 > &ctr, ResultType *r, const std::array< T, K+4 > &par)
Definition: threefry_avx2_2x32.hpp:124
MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_KBOX
#define MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_KBOX(N)
Definition: threefry_avx2_2x32.hpp:42
MCKL_PUSH_GCC_WARNING
#define MCKL_PUSH_GCC_WARNING(warning)
Definition: compiler.h:78
MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_CYCLE_8
#define MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_CYCLE_8(N)
Definition: threefry_avx2_2x32.hpp:89
mckl::internal::Threefry2x32GeneratorAVX2Impl::eval
static void eval(const void *plain, void *cipher, const std::array< T, K+4 > &par)
Definition: threefry_avx2_2x32.hpp:117
mckl::internal::Threefry2x32GeneratorGenericImpl::eval
static void eval(const void *plain, void *cipher, const std::array< T, K+4 > &par)
Definition: threefry_generic_2x32.hpp:98
mckl::internal::Threefry2x32GeneratorAVX2Impl::eval
static void eval(std::array< std::uint64_t, 1 > &ctr, std::size_t n, ResultType *r, const std::array< T, K+4 > &par)
Definition: threefry_avx2_2x32.hpp:131
mckl::ThreefryConstants
Default Threefry constants.
Definition: threefry_constants.hpp:147
threefry_constants.hpp
mckl
Definition: mcmc.hpp:40
threefry_common.hpp
MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_CYCLE_4
#define MCKL_RANDOM_INTERNAL_THREEFRY_AVX2_2X32_CYCLE_4(N)
Definition: threefry_avx2_2x32.hpp:82
MCKL_POP_GCC_WARNING
#define MCKL_POP_GCC_WARNING
Definition: compiler.h:79
Generated by   doxygen 1.8.13