MCKL
Monte Carlo Kernel Library
threefry_sse2_4x32.hpp
Go to the documentation of this file.
1 //============================================================================
2 // MCKL/include/mckl/random/internal/threefry_sse2_4x32.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_SSE2_4X32_HPP
33 #define MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_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_4x32.hpp>
39 
40 MCKL_PUSH_GCC_WARNING("-Wignored-attributes")
41 
42 #define MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_KBOX(N) \
43  xmmt0 = _mm_set1_epi32( \
44  static_cast<int>(ThreefryKBox<T, K, N>::template key<0>(par))); \
45  xmmt1 = _mm_set1_epi32( \
46  static_cast<int>(ThreefryKBox<T, K, N>::template key<1>(par))); \
47  xmmt2 = _mm_set1_epi32( \
48  static_cast<int>(ThreefryKBox<T, K, N>::template key<2>(par))); \
49  xmmt3 = _mm_set1_epi32( \
50  static_cast<int>(ThreefryKBox<T, K, N>::template key<3>(par))); \
51  xmms0 = _mm_add_epi32(xmms0, xmmt0); \
52  xmms1 = _mm_add_epi32(xmms1, xmmt1); \
53  xmms2 = _mm_add_epi32(xmms2, xmmt2); \
54  xmms3 = _mm_add_epi32(xmms3, xmmt3); \
55  xmms4 = _mm_add_epi32(xmms4, xmmt0); \
56  xmms5 = _mm_add_epi32(xmms5, xmmt1); \
57  xmms6 = _mm_add_epi32(xmms6, xmmt2); \
58  xmms7 = _mm_add_epi32(xmms7, xmmt3);
59 
60 #define MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_RBOX(N) \
61  { \
62  constexpr int L0 = Constants::rotate::value[0][(N - 1) % 8]; \
63  constexpr int L1 = Constants::rotate::value[1][(N - 1) % 8]; \
64  constexpr int R0 = 32 - L0; \
65  constexpr int R1 = 32 - L1; \
66  xmms0 = _mm_add_epi32(xmms0, xmms1); \
67  xmms2 = _mm_add_epi32(xmms2, xmms3); \
68  xmms4 = _mm_add_epi32(xmms4, xmms5); \
69  xmms6 = _mm_add_epi32(xmms6, xmms7); \
70  xmmt1 = _mm_slli_epi32(xmms1, L0); \
71  xmmt3 = _mm_slli_epi32(xmms3, L1); \
72  xmmt5 = _mm_slli_epi32(xmms5, L0); \
73  xmmt7 = _mm_slli_epi32(xmms7, L1); \
74  xmms1 = _mm_srli_epi32(xmms1, R0); \
75  xmms3 = _mm_srli_epi32(xmms3, R1); \
76  xmms5 = _mm_srli_epi32(xmms5, R0); \
77  xmms7 = _mm_srli_epi32(xmms7, R1); \
78  xmmt1 = _mm_or_si128(xmms1, xmmt1); \
79  xmmt3 = _mm_or_si128(xmms3, xmmt3); \
80  xmmt5 = _mm_or_si128(xmms5, xmmt5); \
81  xmmt7 = _mm_or_si128(xmms7, xmmt7); \
82  xmms1 = _mm_xor_si128(xmms2, xmmt3); \
83  xmms3 = _mm_xor_si128(xmms0, xmmt1); \
84  xmms5 = _mm_xor_si128(xmms6, xmmt7); \
85  xmms7 = _mm_xor_si128(xmms4, xmmt5); \
86  }
87 
88 #define MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_CYCLE_4(N) \
89  MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_RBOX(N * 8 + 1); \
90  MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_RBOX(N * 8 + 2); \
91  MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_RBOX(N * 8 + 3); \
92  MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_RBOX(N * 8 + 4); \
93  MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_KBOX(N * 8 + 4);
94 
95 #define MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_CYCLE_8(N) \
96  MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_RBOX(N * 8 + 1); \
97  MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_RBOX(N * 8 + 2); \
98  MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_RBOX(N * 8 + 3); \
99  MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_RBOX(N * 8 + 4); \
100  MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_KBOX(N * 8 + 4); \
101  MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_RBOX(N * 8 + 5); \
102  MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_RBOX(N * 8 + 6); \
103  MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_RBOX(N * 8 + 7); \
104  MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_RBOX(N * 8 + 8); \
105  MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_KBOX(N * 8 + 8);
106 
107 namespace mckl {
108 
109 namespace internal {
110 
111 template <typename T>
112 class Threefry4x32GeneratorSSE2Impl
113 {
114  static_assert(std::numeric_limits<T>::digits == 32,
115  "**Threefry4x32GeneratorSSE2Impl** used with T other than a "
116  "32-bit unsigned integers");
117 
118  static constexpr std::size_t K = 4;
119  static constexpr std::size_t Rounds = 20;
120  using Constants = ThreefryConstants<T, 4>;
121 
122  public:
123  static void eval(
124  const void *plain, void *cipher, const std::array<T, K + 4> &par)
125  {
126  Threefry4x32GeneratorGenericImpl<T>::eval(plain, cipher, par);
127  }
128 
129  template <typename ResultType>
130  static void eval(std::array<std::uint64_t, 2> &ctr, ResultType *r,
131  const std::array<T, K + 4> &par)
132  {
133  Threefry4x32GeneratorGenericImpl<T>::eval(ctr, r, par);
134  }
135 
136  template <typename ResultType>
137  static void eval(std::array<std::uint64_t, 2> &ctr, std::size_t n,
138  ResultType *r, const std::array<T, K + 4> &par)
139  {
140  constexpr std::size_t R = sizeof(T) * K / sizeof(ResultType);
141 
142  const std::size_t n0 =
143  static_cast<std::size_t>(std::min(static_cast<std::uint64_t>(n),
144  std::numeric_limits<std::uint64_t>::max() - ctr.front()));
145 
146  eval_kernel(ctr, n0, r, par);
147  n -= n0;
148  r += n0 * R;
149 
150  if (n != 0) {
151  eval(ctr, r, par);
152  n -= 1;
153  r += R;
154  }
155 
156  eval_kernel(ctr, n, r, par);
157  }
158 
159  private:
160  template <typename ResultType>
161  static void eval_kernel(std::array<std::uint64_t, 2> &ctr, std::size_t n,
162  ResultType *r, const std::array<T, K + 4> &par)
163  {
164  constexpr std::size_t S = 8;
165  constexpr std::size_t N = sizeof(__m128i) * S / (sizeof(T) * K);
166 
167  __m128i xmmc =
168  _mm_set_epi64x(static_cast<MCKL_INT64>(std::get<1>(ctr)),
169  static_cast<MCKL_INT64>(std::get<0>(ctr)));
170  ctr.front() += n;
171 
172  __m128i *rptr = reinterpret_cast<__m128i *>(r);
173  while (n != 0) {
174  __m128i xmms0 = _mm_add_epi64(xmmc, _mm_set_epi64x(0, 1));
175  __m128i xmms1 = _mm_add_epi64(xmmc, _mm_set_epi64x(0, 2));
176  __m128i xmms2 = _mm_add_epi64(xmmc, _mm_set_epi64x(0, 3));
177  __m128i xmms3 = _mm_add_epi64(xmmc, _mm_set_epi64x(0, 4));
178  __m128i xmms4 = _mm_add_epi64(xmmc, _mm_set_epi64x(0, 5));
179  __m128i xmms5 = _mm_add_epi64(xmmc, _mm_set_epi64x(0, 6));
180  __m128i xmms6 = _mm_add_epi64(xmmc, _mm_set_epi64x(0, 7));
181  __m128i xmms7 = _mm_add_epi64(xmmc, _mm_set_epi64x(0, 8));
182  xmmc = xmms7;
183 
184  __m128i xmmt0;
185  __m128i xmmt1;
186  __m128i xmmt2;
187  __m128i xmmt3;
188  __m128i xmmt4;
189  __m128i xmmt5;
190  __m128i xmmt6;
191  __m128i xmmt7;
192 
193  xmmt1 = xmms1;
194  xmms1 = xmms2;
195  xmms2 = xmms4;
196  xmms4 = xmmt1;
197 
198  xmmt3 = xmms3;
199  xmms3 = xmms6;
200  xmms6 = xmms5;
201  xmms5 = xmmt3;
202 
203  xmmt0 = _mm_unpacklo_epi32(xmms0, xmms1);
204  xmmt1 = _mm_unpacklo_epi32(xmms2, xmms3);
205  xmmt2 = _mm_unpackhi_epi32(xmms0, xmms1);
206  xmmt3 = _mm_unpackhi_epi32(xmms2, xmms3);
207  xmmt4 = _mm_unpacklo_epi32(xmms4, xmms5);
208  xmmt5 = _mm_unpacklo_epi32(xmms6, xmms7);
209  xmmt6 = _mm_unpackhi_epi32(xmms4, xmms5);
210  xmmt7 = _mm_unpackhi_epi32(xmms6, xmms7);
211  xmms0 = _mm_unpacklo_epi64(xmmt0, xmmt1);
212  xmms1 = _mm_unpackhi_epi64(xmmt0, xmmt1);
213  xmms2 = _mm_unpacklo_epi64(xmmt2, xmmt3);
214  xmms3 = _mm_unpackhi_epi64(xmmt2, xmmt3);
215  xmms4 = _mm_unpacklo_epi64(xmmt4, xmmt5);
216  xmms5 = _mm_unpackhi_epi64(xmmt4, xmmt5);
217  xmms6 = _mm_unpacklo_epi64(xmmt6, xmmt7);
218  xmms7 = _mm_unpackhi_epi64(xmmt6, xmmt7);
219 
220  MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_KBOX(0)
221  MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_CYCLE_8(0)
222  MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_CYCLE_8(1)
223  MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_CYCLE_4(2)
224 
225  xmmt0 = _mm_unpacklo_epi32(xmms0, xmms1);
226  xmmt1 = _mm_unpacklo_epi32(xmms2, xmms3);
227  xmmt2 = _mm_unpackhi_epi32(xmms0, xmms1);
228  xmmt3 = _mm_unpackhi_epi32(xmms2, xmms3);
229  xmmt4 = _mm_unpacklo_epi32(xmms4, xmms5);
230  xmmt5 = _mm_unpacklo_epi32(xmms6, xmms7);
231  xmmt6 = _mm_unpackhi_epi32(xmms4, xmms5);
232  xmmt7 = _mm_unpackhi_epi32(xmms6, xmms7);
233  xmms0 = _mm_unpacklo_epi64(xmmt0, xmmt1);
234  xmms1 = _mm_unpackhi_epi64(xmmt0, xmmt1);
235  xmms2 = _mm_unpacklo_epi64(xmmt2, xmmt3);
236  xmms3 = _mm_unpackhi_epi64(xmmt2, xmmt3);
237  xmms4 = _mm_unpacklo_epi64(xmmt4, xmmt5);
238  xmms5 = _mm_unpackhi_epi64(xmmt4, xmmt5);
239  xmms6 = _mm_unpacklo_epi64(xmmt6, xmmt7);
240  xmms7 = _mm_unpackhi_epi64(xmmt6, xmmt7);
241 
242  xmmt2 = xmms2;
243  xmms2 = xmms1;
244  xmms1 = xmms4;
245  xmms4 = xmmt2;
246 
247  xmmt6 = xmms6;
248  xmms6 = xmms3;
249  xmms3 = xmms5;
250  xmms5 = xmmt6;
251 
252  if (n >= N) {
253  n -= N;
254  _mm_storeu_si128(rptr++, xmms0);
255  _mm_storeu_si128(rptr++, xmms1);
256  _mm_storeu_si128(rptr++, xmms2);
257  _mm_storeu_si128(rptr++, xmms3);
258  _mm_storeu_si128(rptr++, xmms4);
259  _mm_storeu_si128(rptr++, xmms5);
260  _mm_storeu_si128(rptr++, xmms6);
261  _mm_storeu_si128(rptr++, xmms7);
262  } else {
263  std::array<__m128i, S> s;
264  std::get<0>(s) = xmms0;
265  std::get<1>(s) = xmms1;
266  std::get<2>(s) = xmms2;
267  std::get<3>(s) = xmms3;
268  std::get<4>(s) = xmms4;
269  std::get<5>(s) = xmms5;
270  std::get<6>(s) = xmms6;
271  std::get<7>(s) = xmms7;
272  std::memcpy(rptr, s.data(), n * sizeof(T) * K);
273  break;
274  }
275  }
276  }
277 }; // class Threefry4x32GeneratorSSE2Impl
278 
279 } // namespace internal
280 
281 } // namespace mckl
282 
283 MCKL_POP_GCC_WARNING
284 
285 #endif // MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_HPP
threefry_generic_4x32.hpp
MCKL_PUSH_GCC_WARNING
#define MCKL_PUSH_GCC_WARNING(warning)
Definition: compiler.h:78
mckl::internal::Threefry4x32GeneratorSSE2Impl::eval
static void eval(std::array< std::uint64_t, 2 > &ctr, ResultType *r, const std::array< T, K+4 > &par)
Definition: threefry_sse2_4x32.hpp:130
mckl::internal::Threefry4x32GeneratorSSE2Impl::eval
static void eval(const void *plain, void *cipher, const std::array< T, K+4 > &par)
Definition: threefry_sse2_4x32.hpp:123
mckl::internal::Threefry4x32GeneratorGenericImpl::eval
static void eval(const void *plain, void *cipher, const std::array< T, K+4 > &par)
Definition: threefry_generic_4x32.hpp:111
MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_CYCLE_8
#define MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_CYCLE_8(N)
Definition: threefry_sse2_4x32.hpp:95
mckl::ThreefryConstants
Default Threefry constants.
Definition: threefry_constants.hpp:147
MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_KBOX
#define MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_KBOX(N)
Definition: threefry_sse2_4x32.hpp:42
threefry_constants.hpp
mckl
Definition: mcmc.hpp:40
threefry_common.hpp
MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_CYCLE_4
#define MCKL_RANDOM_INTERNAL_THREEFRY_SSE2_4X32_CYCLE_4(N)
Definition: threefry_sse2_4x32.hpp:88
mckl::internal::Threefry4x32GeneratorSSE2Impl::eval
static void eval(std::array< std::uint64_t, 2 > &ctr, std::size_t n, ResultType *r, const std::array< T, K+4 > &par)
Definition: threefry_sse2_4x32.hpp:137
MCKL_POP_GCC_WARNING
#define MCKL_POP_GCC_WARNING
Definition: compiler.h:79
Generated by   doxygen 1.8.13