値の生成方法
生成処理は、Intel IA32 の SSE4.2 に対応したプロセッサ命令で
RDTSC
CRC32 EAX, EDX
二命令だけです。 EAX の値が乱数のようなものになるだろうという期待をします。
RDTSC と CRC32 の関係
CPU の CRC32 は元となる crc32 計算
- eax = crc32(edx)
とは少し違い
- eax = crc32(eax ^ edx)
という計算をします。(レジスタのビット反転を含むので詳しくは仕様書を参照してください)
元の y = CRC32(x) では、複数の x に対して、y が1つ(重複)となることがありません。符号なし32ビット整数における全単射になるでしょうか。(過去記事)
RDTSC は EDX:EAX にプロセッサのタイムスタンプ(分解能は CPU クロック) を 64 ビットで返します。EDX は 232 クロック置きに更新されます。4GHz 動作で約 1 秒間隔です。
これを生成プログラムのように CRC32 に渡すと、RDTSC[EDX:EAX] と CRC32 の関係は
EDX | EAX | CRC32 |
---|---|---|
00000000 |
00000000 |
00000000 |
00000000 |
00000001 |
dd45aab8 |
00000000 |
00000002 |
bf672381 |
00000000 |
00000003 |
62228939 |
... | ... | ... |
00000000 |
fffffffe |
6add1e80 |
00000000 |
ffffffff |
b798b438 |
00000001 |
00000000 |
dd45aab8 |
00000001 |
00000001 |
00000000 |
00000001 |
00000002 |
62228939 |
00000001 |
00000003 |
bf672381 |
... | ... | ... |
fffffffe |
fffffffe |
00000000 |
fffffffe |
ffffffff |
dd45aab8 |
ffffffff |
00000000 |
b798b438 |
ffffffff |
00000001 |
6add1e80 |
... | ... | ... |
ffffffff |
fffffffe |
dd45aab8 |
ffffffff |
ffffffff |
00000000 |
CRC32命令の計算
def crc32(eax, edx):
eax ^= edx
for _ in range(32):
f = eax & 1
eax >>= 1
if not f:
continue
eax ^= 0x82f63b78
return eax
となるので、 TSC = (EDX << 32) | EAX
に対する CRC32[TSC] を乱数列(質に疑問符あり)と捉えます。周期は 264 CPU サイクルで、4GHz 動作だと百数十年になります。
生成プログラムは、呼び出された時点の RDTSC(時間とともに増える)から CRC32[TSC] で引き出すことになります。結果、乱数列臭いものが出てくることは予想できますが、実行状況によっては乱数列とよべないものが並ぶことも考えられます。CPU の動作に詳しい人なら、狙ったタイミングで RDTSC を実行する難しさはご存じでしょう。よって、乱数列とよべないものが並ぶことは滅多にないと考えられますが、確率がゼロでないことは実際に発生しうることを意味します。
Python モジュール spcrand
RDTSC を使う方法なので疑似乱数ではありません。
RDTSC はシステムによって無効になっている場合があります。
どんな数列になるかは実行してみなければ分からないので Python モジュールにしてみました。
メソッドは1つです
spcrand() -> int: 符号なし整数(32ビット)で生成プログラム一回分
spcrand(count: int) -> bytearray: 生成プログラム count 回分で count * 4 バイトの数列
動作確認している環境
OS
macOS 12.2.1
C コンパイラ (Xcode 13.2.1)
Apple clang version 13.0.0 (clang-1300.0.29.30)
Target: x86_64-apple-darwin21.3.0
Python : /usr/bin/python3
Python 3.8.9 (default, Oct 26 2021, 07:25:54)
[Clang 13.0.0 (clang-1300.0.29.30)] on darwin
#ifndef PY_SSIZE_T_CLEAN
#define PY_SSIZE_T_CLEAN
#endif /* PY_SSIZE_T_CLEAN */
#include <Python.h>
#if (defined(i386) || \
defined(__i386) || \
defined(__i386__) || \
defined(__i486__) || \
defined(__i586__) || \
defined(__i686__) || \
defined(__amd64) || \
defined(__amd64__) || \
defined(__x86_64) || \
defined(__x86_64__) || \
0)
#define ISA_X86 1
#else
#define ISA_X86 0
#endif
/*
* Method: spcrand
*/
#if ISA_X86
#define HAS_CHECK_SSE42 1
typedef struct cpuid_result {
unsigned int eax; /* reg:00 */
unsigned int ecx; /* reg:01 */
unsigned int edx; /* reg:10 */
unsigned int ebx; /* reg:11 */
} cpuid_result;
inline static int
check_sse42()
{
cpuid_result res;
res.eax = 1;
res.ecx = 0;
res.edx = 0;
res.ebx = 0;
__asm__ volatile (
"mov %0,%%eax;"
"mov %1,%%ecx;"
"xor %%edx,%%edx;"
"xor %%ebx,%%ebx;"
"cpuid;"
"mov %%eax,%0;"
"mov %%ecx,%1;"
"mov %%edx,%2;"
"mov %%ebx,%3;"
: "+m" (res.eax),
"+m" (res.ecx),
"+m" (res.edx),
"+m" (res.ebx)
: /**/
: "eax", "ecx", "edx", "ebx");
return (res.ecx >> 20) & 1;
}
inline static uint32_t ezrand()
{
#if defined(__GNUC__)
uint32_t eax, edx;
__asm__ volatile (
"rdtsc\n"
"crc32l\t%%edx, %%eax"
: "=a" (eax), "=d" (edx));
return eax;
#else /* !__GNUC__ */
#error "unknown compiler"
#endif
}
#else /* !ISA_X86 */
#error "unknown processor"
#endif
static PyObject *
method_spcrand(PyObject *unused, PyObject *const *args, Py_ssize_t size)
{
PyByteArrayObject *ba;
PyObject *r = NULL;
uint32_t *buffer;
Py_ssize_t count;
Py_ssize_t i;
(void) unused;
if (size == 0)
return PyLong_FromUnsignedLong(ezrand());
if (size != 1) {
PyErr_SetString(PyExc_TypeError, "spcrand expected at 1 argument");
return NULL;
}
if (!PyLong_CheckExact(args[0])) {
PyErr_SetObject(PyExc_TypeError, args[0]);
return NULL;
}
count = PyLong_AsSsize_t(args[0]);
if (count < 0) {
if (!PyErr_Occurred())
PyErr_SetObject(PyExc_TypeError, args[0]);
return NULL;
}
if (!(r = PyByteArray_FromStringAndSize("", 0)))
return NULL;
if (PyByteArray_Resize(r, count * sizeof(uint32_t)) < 0) {
Py_XDECREF(r);
return NULL;
}
ba = (PyByteArrayObject *) r;
buffer = (uint32_t *) ba->ob_start;
for (i = 0; i < count; i++)
buffer[i] = ezrand();
return r;
}
/*
* Module: spcrand
*/
static PyMethodDef spcrand_methods[] = {
{"spcrand", (PyCFunction) method_spcrand, METH_FASTCALL, NULL},
{NULL, NULL, 0, NULL}, /* end */
};
static PyModuleDef spcrand_def = {
PyModuleDef_HEAD_INIT,
.m_name = "spcrand",
.m_doc = NULL,
.m_size = -1,
.m_methods = spcrand_methods,
};
PyMODINIT_FUNC
PyInit_spcrand(void)
{
#if HAS_CHECK_SSE42
if (check_sse42() < 0)
return NULL;
#endif /* HAS_CHECK_SSE42 */
return PyModule_Create(&spcrand_def);
}
#!/usr/bin/env python3
import os
from distutils.core import setup, Extension
def getenv(name, defval=None):
if name in os.environ:
return os.environ[name]
return defval
DEBUG = getenv('DEBUG') in ('true', 'yes')
MAJOR_VERSION = 0
MINOR_VERSION = 1
DEBUG_VERSION = 0
VERSION = '%d.%d.%d' % (MAJOR_VERSION, MINOR_VERSION, DEBUG_VERSION)
DEFINE_MACROS = [
('MAJOR_VERSION', MAJOR_VERSION),
('MINOR_VERSION', MINOR_VERSION),
('DEBUG_VERSION', DEBUG_VERSION),
]
UNDEF_MACROS = []
EXTRA_COMPILE_ARGS = [
'-W',
'-Wall',
'-Wno-invalid-offsetof',
'-Wno-deprecated-declarations',
]
if DEBUG:
DEFINE_MACROS.append(('DEBUG', 1))
UNDEF_MACROS.append('NDEBUG')
EXTRA_COMPILE_ARGS.append('-O0')
pass
setup(name='spcrand',
version=VERSION,
description='',
ext_modules=[Extension(
name='spcrand',
define_macros=DEFINE_MACROS,
undef_macros=UNDEF_MACROS,
extra_compile_args=EXTRA_COMPILE_ARGS,
sources=['spcrand.c'])])
setuptools を使えという警告が出るのはわかっているのですが…
ビルドとテスト
$ env ARCHFLAGS="-arch x86_64" python3 setup.py build
running build
running build_ext
building 'spcrand' extension
clang -Wno-unused-result -Wsign-compare -Wunreachable-code -fno-common -dynamic -DNDEBUG -g -fwrapv -O3 -Wall -iwithsysroot/System/Library/Frameworks/System.framework/PrivateHeaders -iwithsysroot/Applications/Xcode.app/Contents/Developer/Library/Frameworks/Python3.framework/Versions/3.8/Headers -Werror=implicit-function-declaration -arch x86_64 -DMAJOR_VERSION=0 -DMINOR_VERSION=1 -DDEBUG_VERSION=0 -I/Applications/Xcode.app/Contents/Developer/Library/Frameworks/Python3.framework/Versions/3.8/include/python3.8 -c spcrand.c -o build/temp.macosx-10.14-x86_64-3.8/spcrand.o -W -Wall -Wno-invalid-offsetof -Wno-deprecated-declarations
creating build/lib.macosx-10.14-x86_64-3.8
clang -bundle -undefined dynamic_lookup -Wl,-headerpad,0x1000 -arch x86_64 build/temp.macosx-10.14-x86_64-3.8/spcrand.o -o build/lib.macosx-10.14-x86_64-3.8/spcrand.cpython-38-darwin.so
$
$ (cd build/lib.macosx-10.14-x86_64-3.8; python3)
Python 3.8.9 (default, Feb 18 2022, 07:45:34)
[Clang 13.1.6 (clang-1316.0.21.2)] on darwin
Type "help", "copyright", "credits" or "license" for more information.
>>> from spcrand import spcrand
>>> hex(spcrand())
'0x1ff4677b'
>>> hex(spcrand())
'0xe7806714'
>>> hex(spcrand())
'0xf996f420'
>>> r = spcrand(4096)
>>> [r.count(i) for i in range(256)]
[72, 69, 60, 75, 61, 62, 60, 57, 59, 64, 66, 51, 58, 68, 58, 68, 63, 84, 68, 73, 61, 63, 65, 60, 66, 70, 68, 66, 55, 68, 78, 63, 60, 71, 64, 64, 49, 73, 60, 61, 79, 59, 62, 71, 52, 69, 63, 67, 68, 45, 51, 69, 57, 73, 70, 45, 72, 61, 58, 75, 55, 62, 72, 62, 55, 60, 70, 73, 63, 63, 47, 63, 61, 58, 56, 72, 69, 62, 72, 75, 68, 57, 51, 66, 59, 57, 69, 57, 67, 59, 59, 68, 68, 55, 57, 64, 82, 62, 63, 70, 59, 61, 76, 60, 62, 46, 64, 64, 72, 66, 64, 57, 69, 59, 60, 67, 64, 67, 54, 68, 68, 62, 71, 68, 67, 60, 58, 53, 78, 65, 67, 65, 57, 62, 56, 48, 56, 53, 60, 59, 69, 72, 58, 61, 51, 64, 66, 66, 65, 62, 63, 53, 59, 58, 55, 65, 73, 62, 63, 72, 64, 70, 61, 84, 54, 71, 72, 76, 67, 63, 67, 55, 59, 73, 62, 61, 56, 68, 64, 63, 72, 67, 77, 63, 66, 69, 54, 60, 64, 73, 59, 67, 67, 76, 58, 74, 66, 62, 64, 61, 57, 70, 69, 61, 64, 75, 79, 77, 59, 76, 67, 51, 67, 66, 58, 64, 57, 53, 58, 63, 68, 65, 74, 52, 69, 72, 64, 59, 65, 63, 69, 59, 61, 51, 65, 88, 57, 71, 63, 57, 62, 62, 60, 61, 67, 73, 61, 66, 69, 65, 75, 65, 75, 51, 62, 74]
>>> min([r.count(i) for i in range(256)])
45
>>> max([r.count(i) for i in range(256)])
88
>>> len(r)
16384
>>> def dump(l):
... for s in range(0, len(l), 16):
... print(' '.join(['%02x' % d for d in l[s:s+16]]))
...
>>> dump(r[:256])
a2 f8 99 3c a5 4f fc ea 06 b8 db 59 b5 19 3e b1
9d fa 3c d5 59 4d 38 2f 1f 9d 7d 8c a3 0f 7e da
e6 9b 1d 9f 4f 5b 78 44 27 28 bd b3 68 8b 9c 9e
5f 0d ba 1f e6 6c d9 a6 28 ec 9f 34 1f 6a b9 b5
b6 aa dc 6e de d9 19 99 10 59 5f 0b 0f 3c 7b ee
92 3e 1e 52 2e ac 1d 04 8d 5b 3a b7 c2 f8 1b 9a
aa 8b de 6d b5 ee fa 88 09 7c f9 de 66 9d d5 37
4e 7e d7 53 79 f8 f1 d2 36 5b d0 ff f7 e8 70 d3
c0 6e 56 52 8f cd 77 7f f4 ac 56 35 48 3e 55 63
0d aa 36 26 f5 90 b4 68 49 02 b7 3e 35 1f f6 19
1d fc f4 7d 0d 5d f2 1f c3 dd b4 8d f4 5b 92 0c
0c 61 10 42 b0 f3 13 14 4d cd 35 8c e4 0d 50 57
79 0f 35 eb c5 9d 36 bd 5b db 75 e7 78 33 d7 b6
c4 a1 d4 e0 5e df d1 31 90 5f 97 a3 97 08 20 10
09 4e 63 4a b5 dc 60 1c 96 34 c2 4d 2a a6 c1 1b
31 fb a3 75 7e 58 82 58 38 88 c7 fb 0f 0e e1 7a
>>> dump(r[256:256+256])
c1 8e a7 e8 04 05 41 4f b8 97 42 19 c4 8a 03 3e
fc c8 07 38 c1 79 63 d1 84 ed 00 94 38 7f 03 c2
7c d7 82 da 41 66 e6 33 bc 58 c0 ab 7b ab e2 b7
b5 2b a4 25 82 ad 82 a4 aa 4e 80 c0 89 a6 22 91
35 34 21 c7 2a 51 05 22 7c 39 0b a9 c0 ab 08 ff
e2 7f 48 f3 5e ed 4b a5 fc 26 8e 4b c1 97 ea a2
5f d1 a9 f8 fb 5a ee 26 d3 b9 ec 42 e4 3f ca c3
2a bf 8c 51 d2 85 0e 1f e5 03 28 9e aa a0 09 b3
17 f9 2c 81 ab 6b 2f d7 6f dc 2b 2d d3 4e 28 7b
aa 57 cd 8a 1c f2 8c b4 d2 72 ca 26 90 9a c9 0e
5e 1a 8f 9c 69 9c a9 1d 10 85 4c ec 5f 26 6d c1
68 a0 4b 40 c1 60 2e 9b 48 92 af b6 f4 00 ac e0
eb 65 88 05 3b bc 08 2f f5 3c 4e bd d7 e8 0e b1
6b 7a 0d e7 cf f1 4a 39 51 b7 09 63 ed 25 0a 35
4f ee cf db e6 2e aa 00 d1 a8 8c 81 9e 0b ad ac
a6 49 a9 aa 29 65 ca f6 66 c6 eb db b6 1f 6b f1
このときのバイト列の分布としては 16384/256=64, min=64-19, max=64+24 の範囲に収まっているようです。
参考:count 回ループ部アセンブラ出力(抜粋)
for (i = 0; i < count; i++)
buffer[i] = ezrand();
LBB1_14:
rdtsc
crc32l %edx, %eax
movl %eax, (%rcx,%rsi,4)
addq $1, %rsi
cmpq %rsi, %r15
jne LBB1_14
使用例
「地上波アナログ時代の砂嵐風画面」の GIF 動画を作るプログラム
#!/usr/bin/env python3
from mkgif import File as GIF
from spcrand import spcrand
width = 320
height = 240
pixels = width * height
gif = GIF()
gif.append_logical_screen_descriptor(
width, height, global_color=tuple((n, n, n) for n in range(256)))
gif.append_animaition_loop()
for n in range(25):
print('Frame:', n)
gif.append_graphic_control_extension(2)
gif.append_image_descriptr(spcrand() & 0xff for _ in range((pixels)))
pass
gif.append_trailer()
import sys
open(sys.argv[1], 'wb').write(gif)
mkgif は記事「GIF 画像データを作る」にあります。
乱数は利用目的に適合していれば、いわゆる乱数評価方法において不適格なものでも問題ありません。画面上での演出など、簡単な乱数生成で十分な場合が多々あります。
タイトルからして、この記事は「ポンコツ」という意味を含ませています。
よって「この乱数アルゴリズムを使え」みたいなコメントは
その通り過ぎて鬱陶しいだけなので一切不要です。
検証する時は、RDTSC+CRC32 自身に数列を決定する能力がないので、「プログラム・実行環境・実行時」がセットで決定されることに注意しましょう。実行環境には物理的要因(主に機器温度:温度上昇を抑えるためにパフォーマンスを低下させる仕組みが入っている場合)も考えられます。
検証記事は見つけたら追加していきます。