← Back to branch summary

~ppsspp/ppsspp/ppsspp_1.3.0

« back to all changes in this revision

Viewing changes to ext/sfmt19937/SFMT.c

  • Committer: Sérgio Benjamim
  • Date: 2017-01-02 00:12:05 UTC
  • Revision ID: sergio_br2@yahoo.com.br-20170102001205-cxbta9za203nmjwm
1.3.0 source (from ppsspp_1.3.0-r160.p5.l1762.a165.t83~56~ubuntu16.04.1.tar.xz).

Show diffs side-by-side

added added

removed removed

Lines of Context:
ext/sfmt19937/SFMT.c
 
1
/**
 
2
 * @file  SFMT.c
 
3
 * @brief SIMD oriented Fast Mersenne Twister(SFMT)
 
4
 *
 
5
 * @author Mutsuo Saito (Hiroshima University)
 
6
 * @author Makoto Matsumoto (Hiroshima University)
 
7
 *
 
8
 * Copyright (C) 2006, 2007 Mutsuo Saito, Makoto Matsumoto and Hiroshima
 
9
 * University.
 
10
 * Copyright (C) 2012 Mutsuo Saito, Makoto Matsumoto, Hiroshima
 
11
 * University and The University of Tokyo.
 
12
 * Copyright (C) 2013 Mutsuo Saito, Makoto Matsumoto and Hiroshima
 
13
 * University.
 
14
 * All rights reserved.
 
15
 *
 
16
 * The 3-clause BSD License is applied to this software, see
 
17
 * LICENSE.txt
 
18
 */
 
19
 
 
20
#if defined(__cplusplus)
 
21
extern "C" {
 
22
#endif
 
23
 
 
24
#include <string.h>
 
25
#include <assert.h>
 
26
#include "SFMT.h"
 
27
#include "SFMT-params.h"
 
28
#include "SFMT-common.h"
 
29
 
 
30
#if defined(__BIG_ENDIAN__) && !defined(__amd64) && !defined(BIG_ENDIAN64)
 
31
#define BIG_ENDIAN64 1
 
32
#endif
 
33
#if defined(HAVE_ALTIVEC) && !defined(BIG_ENDIAN64)
 
34
#define BIG_ENDIAN64 1
 
35
#endif
 
36
#if defined(ONLY64) && !defined(BIG_ENDIAN64)
 
37
  #if defined(__GNUC__)
 
38
    #error "-DONLY64 must be specified with -DBIG_ENDIAN64"
 
39
  #endif
 
40
#undef ONLY64
 
41
#endif
 
42
 
 
43
/**
 
44
 * parameters used by sse2.
 
45
 */
 
46
static const w128_t sse2_param_mask = {{SFMT_MSK1, SFMT_MSK2,
 
47
                                        SFMT_MSK3, SFMT_MSK4}};
 
48
/*----------------
 
49
  STATIC FUNCTIONS
 
50
  ----------------*/
 
51
inline static int idxof(int i);
 
52
inline static void gen_rand_array(sfmt_t * sfmt, w128_t *array, int size);
 
53
inline static uint32_t func1(uint32_t x);
 
54
inline static uint32_t func2(uint32_t x);
 
55
static void period_certification(sfmt_t * sfmt);
 
56
#if defined(BIG_ENDIAN64) && !defined(ONLY64)
 
57
inline static void swap(w128_t *array, int size);
 
58
#endif
 
59
 
 
60
#if defined(HAVE_ALTIVEC)
 
61
  #include "SFMT-alti.h"
 
62
#elif defined(HAVE_SSE2)
 
63
  #if defined(_MSC_VER)
 
64
    #include "SFMT-sse2-msc.h"
 
65
  #else
 
66
    #include "SFMT-sse2.h"
 
67
  #endif
 
68
#endif
 
69
 
 
70
/**
 
71
 * This function simulate a 64-bit index of LITTLE ENDIAN
 
72
 * in BIG ENDIAN machine.
 
73
 */
 
74
#ifdef ONLY64
 
75
inline static int idxof(int i) {
 
76
    return i ^ 1;
 
77
}
 
78
#else
 
79
inline static int idxof(int i) {
 
80
    return i;
 
81
}
 
82
#endif
 
83
 
 
84
#if (!defined(HAVE_ALTIVEC)) && (!defined(HAVE_SSE2))
 
85
/**
 
86
 * This function fills the user-specified array with pseudorandom
 
87
 * integers.
 
88
 *
 
89
 * @param sfmt SFMT internal state
 
90
 * @param array an 128-bit array to be filled by pseudorandom numbers.
 
91
 * @param size number of 128-bit pseudorandom numbers to be generated.
 
92
 */
 
93
inline static void gen_rand_array(sfmt_t * sfmt, w128_t *array, int size) {
 
94
    int i, j;
 
95
    w128_t *r1, *r2;
 
96
 
 
97
    r1 = &sfmt->state[SFMT_N - 2];
 
98
    r2 = &sfmt->state[SFMT_N - 1];
 
99
    for (i = 0; i < SFMT_N - SFMT_POS1; i++) {
 
100
        do_recursion(&array[i], &sfmt->state[i], &sfmt->state[i + SFMT_POS1], r1, r2);
 
101
        r1 = r2;
 
102
        r2 = &array[i];
 
103
    }
 
104
    for (; i < SFMT_N; i++) {
 
105
        do_recursion(&array[i], &sfmt->state[i],
 
106
                     &array[i + SFMT_POS1 - SFMT_N], r1, r2);
 
107
        r1 = r2;
 
108
        r2 = &array[i];
 
109
    }
 
110
    for (; i < size - SFMT_N; i++) {
 
111
        do_recursion(&array[i], &array[i - SFMT_N],
 
112
                     &array[i + SFMT_POS1 - SFMT_N], r1, r2);
 
113
        r1 = r2;
 
114
        r2 = &array[i];
 
115
    }
 
116
    for (j = 0; j < 2 * SFMT_N - size; j++) {
 
117
        sfmt->state[j] = array[j + size - SFMT_N];
 
118
    }
 
119
    for (; i < size; i++, j++) {
 
120
        do_recursion(&array[i], &array[i - SFMT_N],
 
121
                     &array[i + SFMT_POS1 - SFMT_N], r1, r2);
 
122
        r1 = r2;
 
123
        r2 = &array[i];
 
124
        sfmt->state[j] = array[i];
 
125
    }
 
126
}
 
127
#endif
 
128
 
 
129
#if defined(BIG_ENDIAN64) && !defined(ONLY64) && !defined(HAVE_ALTIVEC)
 
130
inline static void swap(w128_t *array, int size) {
 
131
    int i;
 
132
    uint32_t x, y;
 
133
 
 
134
    for (i = 0; i < size; i++) {
 
135
        x = array[i].u[0];
 
136
        y = array[i].u[2];
 
137
        array[i].u[0] = array[i].u[1];
 
138
        array[i].u[2] = array[i].u[3];
 
139
        array[i].u[1] = x;
 
140
        array[i].u[3] = y;
 
141
    }
 
142
}
 
143
#endif
 
144
/**
 
145
 * This function represents a function used in the initialization
 
146
 * by init_by_array
 
147
 * @param x 32-bit integer
 
148
 * @return 32-bit integer
 
149
 */
 
150
static uint32_t func1(uint32_t x) {
 
151
    return (x ^ (x >> 27)) * (uint32_t)1664525UL;
 
152
}
 
153
 
 
154
/**
 
155
 * This function represents a function used in the initialization
 
156
 * by init_by_array
 
157
 * @param x 32-bit integer
 
158
 * @return 32-bit integer
 
159
 */
 
160
static uint32_t func2(uint32_t x) {
 
161
    return (x ^ (x >> 27)) * (uint32_t)1566083941UL;
 
162
}
 
163
 
 
164
/**
 
165
 * This function certificate the period of 2^{MEXP}
 
166
 * @param sfmt SFMT internal state
 
167
 */
 
168
static void period_certification(sfmt_t * sfmt) {
 
169
    int inner = 0;
 
170
    int i, j;
 
171
    uint32_t work;
 
172
    uint32_t *psfmt32 = &sfmt->state[0].u[0];
 
173
    const uint32_t parity[4] = {SFMT_PARITY1, SFMT_PARITY2,
 
174
                                SFMT_PARITY3, SFMT_PARITY4};
 
175
 
 
176
    for (i = 0; i < 4; i++)
 
177
        inner ^= psfmt32[idxof(i)] & parity[i];
 
178
    for (i = 16; i > 0; i >>= 1)
 
179
        inner ^= inner >> i;
 
180
    inner &= 1;
 
181
    /* check OK */
 
182
    if (inner == 1) {
 
183
        return;
 
184
    }
 
185
    /* check NG, and modification */
 
186
    for (i = 0; i < 4; i++) {
 
187
        work = 1;
 
188
        for (j = 0; j < 32; j++) {
 
189
            if ((work & parity[i]) != 0) {
 
190
                psfmt32[idxof(i)] ^= work;
 
191
                return;
 
192
            }
 
193
            work = work << 1;
 
194
        }
 
195
    }
 
196
}
 
197
 
 
198
/*----------------
 
199
  PUBLIC FUNCTIONS
 
200
  ----------------*/
 
201
#define UNUSED_VARIABLE(x) (void)(x)
 
202
/**
 
203
 * This function returns the identification string.
 
204
 * The string shows the word size, the Mersenne exponent,
 
205
 * and all parameters of this generator.
 
206
 * @param sfmt SFMT internal state
 
207
 */
 
208
const char *sfmt_get_idstring(sfmt_t * sfmt) {
 
209
    UNUSED_VARIABLE(sfmt);
 
210
    return SFMT_IDSTR;
 
211
}
 
212
 
 
213
/**
 
214
 * This function returns the minimum size of array used for \b
 
215
 * fill_array32() function.
 
216
 * @param sfmt SFMT internal state
 
217
 * @return minimum size of array used for fill_array32() function.
 
218
 */
 
219
int sfmt_get_min_array_size32(sfmt_t * sfmt) {
 
220
    UNUSED_VARIABLE(sfmt);
 
221
    return SFMT_N32;
 
222
}
 
223
 
 
224
/**
 
225
 * This function returns the minimum size of array used for \b
 
226
 * fill_array64() function.
 
227
 * @param sfmt SFMT internal state
 
228
 * @return minimum size of array used for fill_array64() function.
 
229
 */
 
230
int sfmt_get_min_array_size64(sfmt_t * sfmt) {
 
231
    UNUSED_VARIABLE(sfmt);
 
232
    return SFMT_N64;
 
233
}
 
234
 
 
235
#if !defined(HAVE_SSE2) && !defined(HAVE_ALTIVEC)
 
236
/**
 
237
 * This function fills the internal state array with pseudorandom
 
238
 * integers.
 
239
 * @param sfmt SFMT internal state
 
240
 */
 
241
void sfmt_gen_rand_all(sfmt_t * sfmt) {
 
242
    int i;
 
243
    w128_t *r1, *r2;
 
244
 
 
245
    r1 = &sfmt->state[SFMT_N - 2];
 
246
    r2 = &sfmt->state[SFMT_N - 1];
 
247
    for (i = 0; i < SFMT_N - SFMT_POS1; i++) {
 
248
        do_recursion(&sfmt->state[i], &sfmt->state[i],
 
249
                     &sfmt->state[i + SFMT_POS1], r1, r2);
 
250
        r1 = r2;
 
251
        r2 = &sfmt->state[i];
 
252
    }
 
253
    for (; i < SFMT_N; i++) {
 
254
        do_recursion(&sfmt->state[i], &sfmt->state[i],
 
255
                     &sfmt->state[i + SFMT_POS1 - SFMT_N], r1, r2);
 
256
        r1 = r2;
 
257
        r2 = &sfmt->state[i];
 
258
    }
 
259
}
 
260
#endif
 
261
 
 
262
#ifndef ONLY64
 
263
/**
 
264
 * This function generates pseudorandom 32-bit integers in the
 
265
 * specified array[] by one call. The number of pseudorandom integers
 
266
 * is specified by the argument size, which must be at least 624 and a
 
267
 * multiple of four.  The generation by this function is much faster
 
268
 * than the following gen_rand function.
 
269
 *
 
270
 * For initialization, init_gen_rand or init_by_array must be called
 
271
 * before the first call of this function. This function can not be
 
272
 * used after calling gen_rand function, without initialization.
 
273
 *
 
274
 * @param sfmt SFMT internal state
 
275
 * @param array an array where pseudorandom 32-bit integers are filled
 
276
 * by this function.  The pointer to the array must be \b "aligned"
 
277
 * (namely, must be a multiple of 16) in the SIMD version, since it
 
278
 * refers to the address of a 128-bit integer.  In the standard C
 
279
 * version, the pointer is arbitrary.
 
280
 *
 
281
 * @param size the number of 32-bit pseudorandom integers to be
 
282
 * generated.  size must be a multiple of 4, and greater than or equal
 
283
 * to (MEXP / 128 + 1) * 4.
 
284
 *
 
285
 * @note \b memalign or \b posix_memalign is available to get aligned
 
286
 * memory. Mac OSX doesn't have these functions, but \b malloc of OSX
 
287
 * returns the pointer to the aligned memory block.
 
288
 */
 
289
void sfmt_fill_array32(sfmt_t * sfmt, uint32_t *array, int size) {
 
290
    assert(sfmt->idx == SFMT_N32);
 
291
    assert(size % 4 == 0);
 
292
    assert(size >= SFMT_N32);
 
293
 
 
294
    gen_rand_array(sfmt, (w128_t *)array, size / 4);
 
295
    sfmt->idx = SFMT_N32;
 
296
}
 
297
#endif
 
298
 
 
299
/**
 
300
 * This function generates pseudorandom 64-bit integers in the
 
301
 * specified array[] by one call. The number of pseudorandom integers
 
302
 * is specified by the argument size, which must be at least 312 and a
 
303
 * multiple of two.  The generation by this function is much faster
 
304
 * than the following gen_rand function.
 
305
 *
 
306
 * @param sfmt SFMT internal state
 
307
 * For initialization, init_gen_rand or init_by_array must be called
 
308
 * before the first call of this function. This function can not be
 
309
 * used after calling gen_rand function, without initialization.
 
310
 *
 
311
 * @param array an array where pseudorandom 64-bit integers are filled
 
312
 * by this function.  The pointer to the array must be "aligned"
 
313
 * (namely, must be a multiple of 16) in the SIMD version, since it
 
314
 * refers to the address of a 128-bit integer.  In the standard C
 
315
 * version, the pointer is arbitrary.
 
316
 *
 
317
 * @param size the number of 64-bit pseudorandom integers to be
 
318
 * generated.  size must be a multiple of 2, and greater than or equal
 
319
 * to (MEXP / 128 + 1) * 2
 
320
 *
 
321
 * @note \b memalign or \b posix_memalign is available to get aligned
 
322
 * memory. Mac OSX doesn't have these functions, but \b malloc of OSX
 
323
 * returns the pointer to the aligned memory block.
 
324
 */
 
325
void sfmt_fill_array64(sfmt_t * sfmt, uint64_t *array, int size) {
 
326
    assert(sfmt->idx == SFMT_N32);
 
327
    assert(size % 2 == 0);
 
328
    assert(size >= SFMT_N64);
 
329
 
 
330
    gen_rand_array(sfmt, (w128_t *)array, size / 2);
 
331
    sfmt->idx = SFMT_N32;
 
332
 
 
333
#if defined(BIG_ENDIAN64) && !defined(ONLY64)
 
334
    swap((w128_t *)array, size /2);
 
335
#endif
 
336
}
 
337
 
 
338
/**
 
339
 * This function initializes the internal state array with a 32-bit
 
340
 * integer seed.
 
341
 *
 
342
 * @param sfmt SFMT internal state
 
343
 * @param seed a 32-bit integer used as the seed.
 
344
 */
 
345
void sfmt_init_gen_rand(sfmt_t * sfmt, uint32_t seed) {
 
346
    int i;
 
347
 
 
348
    uint32_t *psfmt32 = &sfmt->state[0].u[0];
 
349
 
 
350
    psfmt32[idxof(0)] = seed;
 
351
    for (i = 1; i < SFMT_N32; i++) {
 
352
        psfmt32[idxof(i)] = 1812433253UL * (psfmt32[idxof(i - 1)]
 
353
                                            ^ (psfmt32[idxof(i - 1)] >> 30))
 
354
            + i;
 
355
    }
 
356
    sfmt->idx = SFMT_N32;
 
357
    period_certification(sfmt);
 
358
}
 
359
 
 
360
/**
 
361
 * This function initializes the internal state array,
 
362
 * with an array of 32-bit integers used as the seeds
 
363
 * @param sfmt SFMT internal state
 
364
 * @param init_key the array of 32-bit integers, used as a seed.
 
365
 * @param key_length the length of init_key.
 
366
 */
 
367
void sfmt_init_by_array(sfmt_t * sfmt, uint32_t *init_key, int key_length) {
 
368
    int i, j, count;
 
369
    uint32_t r;
 
370
    int lag;
 
371
    int mid;
 
372
    int size = SFMT_N * 4;
 
373
    uint32_t *psfmt32 = &sfmt->state[0].u[0];
 
374
 
 
375
    if (size >= 623) {
 
376
        lag = 11;
 
377
    } else if (size >= 68) {
 
378
        lag = 7;
 
379
    } else if (size >= 39) {
 
380
        lag = 5;
 
381
    } else {
 
382
        lag = 3;
 
383
    }
 
384
    mid = (size - lag) / 2;
 
385
 
 
386
    memset(sfmt, 0x8b, sizeof(sfmt_t));
 
387
    if (key_length + 1 > SFMT_N32) {
 
388
        count = key_length + 1;
 
389
    } else {
 
390
        count = SFMT_N32;
 
391
    }
 
392
    r = func1(psfmt32[idxof(0)] ^ psfmt32[idxof(mid)]
 
393
              ^ psfmt32[idxof(SFMT_N32 - 1)]);
 
394
    psfmt32[idxof(mid)] += r;
 
395
    r += key_length;
 
396
    psfmt32[idxof(mid + lag)] += r;
 
397
    psfmt32[idxof(0)] = r;
 
398
 
 
399
    count--;
 
400
    for (i = 1, j = 0; (j < count) && (j < key_length); j++) {
 
401
        r = func1(psfmt32[idxof(i)] ^ psfmt32[idxof((i + mid) % SFMT_N32)]
 
402
                  ^ psfmt32[idxof((i + SFMT_N32 - 1) % SFMT_N32)]);
 
403
        psfmt32[idxof((i + mid) % SFMT_N32)] += r;
 
404
        r += init_key[j] + i;
 
405
        psfmt32[idxof((i + mid + lag) % SFMT_N32)] += r;
 
406
        psfmt32[idxof(i)] = r;
 
407
        i = (i + 1) % SFMT_N32;
 
408
    }
 
409
    for (; j < count; j++) {
 
410
        r = func1(psfmt32[idxof(i)] ^ psfmt32[idxof((i + mid) % SFMT_N32)]
 
411
                  ^ psfmt32[idxof((i + SFMT_N32 - 1) % SFMT_N32)]);
 
412
        psfmt32[idxof((i + mid) % SFMT_N32)] += r;
 
413
        r += i;
 
414
        psfmt32[idxof((i + mid + lag) % SFMT_N32)] += r;
 
415
        psfmt32[idxof(i)] = r;
 
416
        i = (i + 1) % SFMT_N32;
 
417
    }
 
418
    for (j = 0; j < SFMT_N32; j++) {
 
419
        r = func2(psfmt32[idxof(i)] + psfmt32[idxof((i + mid) % SFMT_N32)]
 
420
                  + psfmt32[idxof((i + SFMT_N32 - 1) % SFMT_N32)]);
 
421
        psfmt32[idxof((i + mid) % SFMT_N32)] ^= r;
 
422
        r -= i;
 
423
        psfmt32[idxof((i + mid + lag) % SFMT_N32)] ^= r;
 
424
        psfmt32[idxof(i)] = r;
 
425
        i = (i + 1) % SFMT_N32;
 
426
    }
 
427
 
 
428
    sfmt->idx = SFMT_N32;
 
429
    period_certification(sfmt);
 
430
}
 
431
#if defined(__cplusplus)
 
432
}
 
433
#endif