IchigoJam の乱数のひみつ

# はじめに

[IchigoJam](https://ichigojam.net/)には、乱数を生成する`RND(n)`関数がある。
また、最近のバージョンでは、乱数系列を初期化する`SRND n`命令が使える。

[IchigoJam BASIC 1.4 コマンド一覧](https://fukuno.jig.jp/app/csv/ichigojam-cmd.html) を見ると、

> RND [計算]
> `PRINT RND(6)`
> 書式：RND(数)
> 優先順位：1
> 0から数未満のでたらめな整数を返す

<!-- separator -->

> SRND [プログラム]
> `SRND 0`
> 書式：SRND 数
> 種を指定して乱数を初期化する

となっており、詳細な仕様は開示されていない。

また、FU-SENさんの「IchigoJam BASIC コマンド一覧」を見ても、詳細な仕様はわからない。
[IchigoJam-BASIC/RND.txt at master · fu-sen/IchigoJam-BASIC · GitHub](https://github.com/fu-sen/IchigoJam-BASIC/blob/master/RND.txt)
[IchigoJam-BASIC/SRND.txt at master · fu-sen/IchigoJam-BASIC · GitHub](https://github.com/fu-sen/IchigoJam-BASIC/blob/master/SRND.txt)

そこで、今回はこの乱数の仕様を調査していく。

# SRNDの効力は？

`SRND`を実行しても、一旦インタラクティブを挟むと乱数の系列が変わってしまうようである。
`SRND`で系列を固定するには、インタラクティブに帰らないように実行させないといけないようだ。

実行例

```
IchigoJam BASIC 1.4.1 by jig.jp
OK
SRND0
OK
?RND(100)
6
OK
SRND0
OK
?RND(100)
58
OK
SRND0:?RND(100)
13
OK
SRND0:?RND(100)
13
OK
```


# 指定の数未満にする方法は？

典型的な擬似乱数は、0以上ライブラリで固定された値(例えばC言語の`RAND_MAX`)未満の整数が返る。
この値から0以上指定の数未満の乱数を作るには、典型的には以下の方法がある。

* 剰余をとる `rand() % 100`
* 縮小する `(int)(rand() * 100.0 / RAND_MAX)`

もし「縮小する」が使われているなら、「指定の数」を少し変えた時、結果も少ししか変わらないはずである。
実験をしてみよう。

```
IchigoJam BASIC 1.4.1 by jig.jp
OK
SRND0:?RND(1000)
613
OK
SRND0:?RND(1001)
320
OK
```

「指定の数」を1しか変えていないのに、結果が約半分になってしまった。
したがって、「縮小する」ではなさそうである。

次に、16進数を用いて、剰余になっていそうかを試してみる。

```
IchigoJam BASIC 1.4.1 by jig.jp
OK
SRND0:?HEX$(RND(#10))
D
OK
SRND0:?HEX$(RND(#100))
CD
OK
SRND0:?HEX$(RND(#1000))
9CD
OK
```

剰余を使っている可能性がありそうだ。
とりあえず、剰余を使っていると仮定して先に進もう。
また、
[いつからその方法で偏りのない乱数が得られると錯覚していた？ - アスペ日記](https://takeda25.hatenablog.jp/entry/20130818/1376837093)
で解説されているように、「中途半端な値が出たらやり直す」という処理が入っているかもしれない。
これも後で考えることにしよう。

# 乱数はバージョンによって違う？

以下の、0～19の各乱数系列についてそれぞれ乱数を5個出力するプログラムを試してみた。

```
10 ' ﾗﾝｽｳ ｲﾁﾗﾝ
20 FOR I=0 TO 19
30 SRND I
40 ?DEC$(I,2);
50 FOR J=0 TO 4
60 ?" ";HEX$(RND(#4000),4);
70 NEXT
80 ?""
90 NEXT
```

その結果、今回試した以下のバージョン全てにおいて、以下の結果が得られた。

* IchigoJam S (1.3.1)
* IchigoKamuy (1.4.1)
* IchigoJam R (1.5b)
* IchigoJam web (1.4.3web)

```
 0 09CD 2DD5 0FED 26FE 2DC4
 1 0DC9 29D1 0BE9 22FA 29E4
 2 01C4 25DC 07E4 2EF7 2584
 3 05C0 21D8 03E0 2AF3 21A4
 4 19DF 3DC7 1FFF 36EC 3D44
 5 1DDB 39C3 1BFB 32E8 3964
 6 11D6 35CE 17F6 3EE5 3504
 7 15D2 31CA 13F2 3AE1 3124
 8 29E9 0DF1 2FC9 06DA 0CC4
 9 2DED 09F5 2BCD 02DE 08E4
10 21E0 05F8 27C0 0ED3 0484
11 25E4 01FC 23C4 0AD7 00A4
12 39FB 1DE3 3FDB 16C8 1C44
13 3DFF 19E7 3BDF 12CC 1864
14 31F2 15EA 37D2 1EC1 1404
15 35F6 11EE 33D6 1AC5 1024
16 0985 2D9D 0FA5 26B6 2FC4
17 0D81 2999 0BA1 22B2 2BE4
18 018C 2594 07AC 2EBF 2784
19 0588 2190 03A8 2ABB 23A4
```

バージョンによらず同じ乱数が得られると予想できる。
さらに、よく見るとなんかパターンがわかりそうなようなわからないような…？

バージョンにはよらなそうなので、今後は IchigoKamuy (1.4.1) で解析を勧めていく。

# 剰余をとるもとの値は？

もとの値を指定の値で割った余りをとることで乱数を得ていると仮定したとき、もとの値を知りたい。

やり直しが入る可能性があるものの、乱数系列と何回目の乱数かを固定すれば、
常に同じ値を指定の値で割った値が得られそうだと考えられる。
同じ値を複数の値で割った余りが得られると、中国剰余定理を用いてもとの値が復元できるはずである。
中国剰余定理の計算は、[以前実装したもの](https://qiita.com/mikecat_mixc/items/b1143557abb05d0f0f3c#cant-restore-the-flag)が使える。

まずはいくつかの値を指定し、乱数を得る。

```
IchigoJam BASIC 1.4.1 by jig.jp
OK
SRND0:?RND(32653)
27492
OK
SRND0:?RND(32687)
14041
OK
SRND0:?RND(32693)
5905
OK
SRND0:?RND(32707)
19628
OK
SRND0:?RND(32713)
11498
OK
SRND0:?RND(32717)
6078
OK
SRND0:?RND(32719)
3368
OK
SRND0:?RND(32749)
28216
OK
```

これらを(剰余,割る値)のペアのリストにまとめると、

```
[(27492,32653),(14041,32687),(5905,32693),(19628,32707),(11498,32713),(6078,32717),(3368,32719),(28216,32749)]
```

となる。これを中国剰余定理で計算すると、もとの値は 44337613 (0x2a489cd) と求まった。

`SRND 0`の系列の2番目の乱数も、プログラムを書いて求めてみる。

```
10 LET[0],32653,32687,32693,32707,32713,32717,32719,32749
20 FOR I=0 TO 7
30 SRND 0:X=RND(9)
40 ?"(";RND([I]);",";[I];"),";
50 NEXT
```

実行結果

```
(30323,32653),(29163,32687),(23556,32693),(10893,32707),(5646,32713),(2208,32717),(507,32719),(9181,32749),
```

もとの値は 1633529301 (0x615dadd5) と求まった。

# 状態空間は？格納場所は？

いくつかの剰余をとる前の乱数の値が求まった。
乱数は32bitっぽい感じだが、乱数が32bitだとしても状態空間が32bitだとは限らない。

メモリ上に乱数の値が格納されている可能性があるので、探してみよう。
乱数の値がメモリ上に見つかれば、まわりに状態空間の残りの情報もあるかもしれない。

[IchigoJam BASIC RAM 4kbyte のつかいかた（メモリマップ） #KidsIT #IchigoJam #lpc1114 / 福野泰介の一日一創 / Create every day by Taisuke Fukuno](https://fukuno.jig.jp/2386)
を参考に、0x10000000 を起点として指定のオフセットを加えた場所から256バイトのデータを
仮想アドレス`# 700`にコピーし、出力するプログラムを用意した。
取得前に乱数系列を設定して乱数を1個取得することで、メモリに乱数が乗ることを期待する。

```
10 ' RAM ｼｭﾄｸ
20 LET[0],#2301,#071B,#181B,#2207,#0212,#1889,#2201,#0212,#7818,#7008,#3301,#3101,#3A01,#D1F9,#4770
30 SRND 0:X=RND(100)
40 INPUT "ｵﾌｾｯﾄ?",A
50 X=USR(#800,A)
60 FOR I=0 TO #FF
70 ?HEX$(PEEK(#700+I),2);
80 IF I%16<>15 ?" "; ELSE ?""
90 NEXT
```

[IchigoJam 1.3β7で新登場「Segmentation Fault」、例外対応で広がった拡張可能性！ #asm #IchigoJam #nepal / 福野泰介の一日一創 / Create every day by Taisuke Fukuno](https://fukuno.jig.jp/2487)
を参考にすると、最近のバージョンではシステムのデータはRAMの後ろの方にあると推測できる。

例えば、オフセット `# F00`を指定し、RAMの最後の部分を見てみた。

```
2C 0A 00 10 00 00 07 02 01 00 00 00 00 6A 00 00
00 00 00 00 12 01 7B EE 00 17 00 00 01 00 00 01
00 00 00 00 00 01 02 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 20 00 00 00 18 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
46 0A 00 00 00 00 00 00 00 00 00 00 00 6C DC 02
00 00 00 12 22 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 E5 55 9A 15 B5 3B 12 1F
33 13 49 05 9A 13 49 05 E0 06 00 10 83 0D 00 10
3D 0A 00 10 3E 0A 00 10 80 06 00 10 3E 0A 00 10
73 E5 3E 00 00 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
```

うーん…先ほど確認した乱数のもとの値はいずれも入っていないようである。
RAMの他の部分も見てみたが、乱数のもとの値は発見できなかった。

30行目を

```
30 SRND 0:X=ABS(100)
```

とし、乱数を取得しないようにしてみた。
`ABS`を入れているのは文字数を合わせて条件をなるべく変えないようにするためである。

同様にオフセット`# F00`を指定した結果は

```
2C 0A 00 10 00 00 07 02 01 00 00 00 00 6A 00 00
00 00 00 00 12 01 D8 9F 00 17 00 00 01 00 00 01
00 00 00 00 00 01 02 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 20 00 00 00 18 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
46 0A 00 00 00 00 00 00 00 00 00 00 00 6C DC 02
00 00 00 12 22 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 E5 55 9A 15
B5 3B 12 1F 33 13 49 05 80 08 00 10 83 0D 00 10
3D 0A 00 10 3E 0A 00 10 80 06 00 10 3E 0A 00 10
71 3C A5 00 00 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
```

ん？さっきと比べて`E5 55 9A 15 B5 3B 12 1F 33 13 49 05`というパターンの位置が4バイトずれている。

乱数系列を変えてみよう。30行目を

```
30 SRND 1:X=ABS(100)
```

とし、同様にオフセット`# F00`を指定する。結果は

```
2C 0A 00 10 00 00 07 02 01 00 00 00 00 6A 00 00
00 00 00 00 12 01 B5 C6 00 17 00 00 01 00 00 01
00 00 00 00 00 01 02 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 20 00 00 00 18 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
46 0A 00 00 00 00 00 00 00 00 00 00 00 6C DC 02
00 00 00 12 22 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 01 00 00 00 E5 55 9A 15
B5 3B 12 1F 33 13 49 05 60 07 00 10 83 0D 00 10
3D 0A 00 10 3E 0A 00 10 80 06 00 10 3E 0A 00 10
AC 6B 5B 00 00 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
```

`E5 55 9A 15 B5 3B 12 1F 33 13 49 05`の前のデータが`00 00 00 00`から`01 00 00 00`になった。
ここに`SRND`で指定した値が入るのかもしれない。

30行目を

```
30 SRND #BEEF
```

とし、同様にオフセット`# F00`を指定すると、結果は

```
2C 0A 00 10 00 00 07 02 01 00 00 00 00 6A 00 00
00 00 00 00 0A 01 30 77 00 17 00 00 01 00 00 01
00 00 00 00 00 01 02 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 20 00 00 00 18 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
46 0A 00 00 00 00 00 00 00 00 00 00 00 6C DC 02
00 00 00 12 22 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 EF BE FF FF E5 55 9A 15
B5 3B 12 1F 33 13 49 05 60 08 00 10 83 0D 00 10
35 0A 00 10 36 0A 00 10 80 06 00 10 36 0A 00 10
B8 BF 9A 00 00 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00
```

となった。どうやらこの場所に`SRND`で指定した値が符号拡張されて入るようだ。

# この値って？

さて、乱数を1個生成すると4バイトずれるようであるパターン
`E5 55 9A 15 B5 3B 12 1F 33 13 49 05`を4バイトずつ10進数にしてみると、
`0x159a55e5 = 362436069, 0x1f123bb5 = 521288629, 0x05491333 = 88675123`となる。
乱数関係で、この値といえば…
[Xorshift - Wikipedia](https://ja.wikipedia.org/w/index.php?title=Xorshift&oldid=79799391)
に載っている`xor128`の`y, z, w`の値にそれぞれなっている！
ということは、`SRND`は`xor128`の`x`を指定した値にし、
`y, z, w`をこれらの値に設定するという動作になっていると予想できる。

本当にXorshiftが使われているのか、試してみよう。
まず、適当に乱数を出力するプログラムを用意する。

```
10 ' ﾗﾝｽｳ ﾃｽﾄ
20 INPUT "ｼｰﾄﾞ?",S
30 SRND S
40 FOR I=0 TO 19
50 FOR J=0 TO 4
60 ?DEC$(RND(#7FFF),6);
70 NEXT
80 ?""
90 NEXT
```

シードに`0`を指定した時の出力は

```
  3862 28817 23783 10341 31427
 28768 17208  5212 30730 26758
 10396 28723  6606 28846 10760
 15625 28866  2101  4861 28223
 26797 10169   124  6045 18667
 13713  6134 29440 14295 13830
 22496  2487 19183 13317  8598
 24442 14793 20827  7463 22474
 22358 22888 11557 23760 14514
 12503 25973 19548 17653  8101
 12123 10811 11018 23906 32612
 17319  9595 24906  2513 29018
 10980  4693  9686 27995  9429
 21814  3501 22153 20335 15579
 31676 17659 18705 16419 32000
  1482 16086 11524 28736  4282
 14017   409 18955  6993 20631
 31604 17796 29544 10792 22614
 10946  3863 20139  3208  7284
 24419 19363    51 14622 15347
```

手元にあったXorshiftの実装を用い、C言語で同様(仮)のプログラムを書く。

```c
#include <stdio.h>

static unsigned int x=123456789;
static unsigned int y=362436069;
static unsigned int z=521288629;
static unsigned int w=88675123;

unsigned int randint(void) {
	unsigned int t;
	t=(x^(x<<11));
	x=y;y=z;z=w;
	w=(w^(w>>19))^(t^(t>>8));
	return w;
}

int main(void) {
	int s, i, j;
	if (scanf("%i", &s) != 1) return 1;
	x = s;
	for (i = 0; i < 20; i++) {
		for (j = 0; j < 5; j++) printf("%6u", randint() % 0x7fff);
		putchar('\n');
	}
	return 0;
}
```

シードに`0`を指定した時の実行結果は

```
  7724 24868 14799 20683 30087
 24770  1649 10425 28694 20749
 20793 24679 13212 24925 21521
 31251 24965  4203  9723 23680
 20827 20338   249 12091  4567
 27427 12268 26114 28590 27660
 12225  4974  5599 26634 17196
 16117 29586  8888 14926 12182
 11949 13009 23115 14754 29028
 25007 19179  6329  2539 16202
 24247 21623 22036 15045 32458
  1871 19191 17046  5027 25270
 21961  9386 19372 23224 18859
 10862  7003 11539  7903 31158
 30586  2552  4643    72 31234
  2964 32172 23049 24705  8565
 28034   819  5144 13987  8495
 30442  2826 26322 21585 12462
 21892  7727  7511  6417 14569
 16071  5960   103 29244 30694
```

残念ながら一致しない！
でも諦めるのはまだ早い。内部状態をXorshiftと比べてみよう。
以下のプログラムで、初期化直後および乱数を生成した後の内部状態を出力する。

```
10 ' ﾗﾝｽｳ ﾅｲﾌﾞ ｼﾞｮｳﾀｲ ｼｭﾄｸ
20 POKE#700,#08,#23,#1B,#02,#5B,#18,#1A,#68,#04,#33,#10,#68,#18,#60,#50,#68,#58,#60,#90,#68,#98,#60,#D0,#68,#D8,#60,#70,#47
30 LET[0],#0FA8,#1000
40 SRND #1234
50 FOR I=0 TO 20
60 X=USR(#700,0)
70 ?HEX$([3],4);HEX$([2],4);" ";HEX$([5],4);HEX$([4],4);" ";HEX$([7],4);HEX$([6],4);" ";HEX$([9],4);HEX$([8],4)
80 X=RND(100)
90 NEXT
```

実行結果は

```
00001234 159A55E5 1F123BB5 05491333
159A55E5 1F123BB5 05491333 05D8301C
1F123BB5 05491333 05D8301C C22A783F
05491333 05D8301C C22A783F 4C6B3C5C
05D8301C C22A783F 4C6B3C5C 01F66E69
C22A783F 4C6B3C5C 01F66E69 C56AE69B
4C6B3C5C 01F66E69 C56AE69B 54109589
01F66E69 C56AE69B 54109589 418CCA8B
C56AE69B 54109589 418CCA8B F3BB61F5
54109589 418CCA8B F3BB61F5 61771F27
418CCA8B F3BB61F5 61771F27 B11B725D
F3BB61F5 61771F27 B11B725D 96E42E67
61771F27 B11B725D 96E42E67 BE784187
B11B725D 96E42E67 BE784187 672FFF48
96E42E67 BE784187 672FFF48 0DCCE06A
BE784187 672FFF48 0DCCE06A BAEC60A2
672FFF48 0DCCE06A BAEC60A2 C6E47A01
0DCCE06A BAEC60A2 C6E47A01 DE29082A
BAEC60A2 C6E47A01 DE29082A B48C6C35
C6E47A01 DE29082A B48C6C35 6DBCE376
DE29082A B48C6C35 6DBCE376 886DA8B2
```

C言語版は

```c
#include <stdio.h>

static unsigned int x=123456789;
static unsigned int y=362436069;
static unsigned int z=521288629;
static unsigned int w=88675123;

unsigned int randint(void) {
	unsigned int t;
	t=(x^(x<<11));
	x=y;y=z;z=w;
	w=(w^(w>>19))^(t^(t>>8));
	return w;
}

int main(void) {
	int i;
	x = 0x1234;
	for (i = 0; i <= 20; i++) {
		printf("%08X %08X %08X %08X\n", x, y, z, w);
		randint();
	}
	return 0;
}
```

このプログラムの実行結果は、一致した！
やはりXorshiftが使われており、1回の`RND(n)`の処理で1回Xorshiftの乱数が使われることが多そうだ。
となると、「Xorshiftの結果」から「乱数のもとの値(剰余をとる時の割られる数)」を求める方法の調査が求められる。

# 乱数のもとの値の求め方？

まず、Xorshiftの計算方法をよく見てみる。

```c
unsigned int randint(void) {
	unsigned int t;
	t=(x^(x<<11));
	x=y;y=z;z=w;
	w=(w^(w>>19))^(t^(t>>8));
	return w;
}
```

すると、`x`に`0`を、`w`に19ビット以下の値を設定しておけば、設定した`w`がそのまま`randint()`の返り値になることがわかる。
さらに、`y`と`z`も`0`に設定しておけば、`randint()`は3回連続で同じ値を返す。
これは中国剰余定理でもとの値を求めるために、同じ値を3種類の値で割った余りを求めるのに都合が良い。

実験をしてみる。

```
10 ' ﾗﾝｽｳ ﾅｲﾌﾞ ｼﾞｮｳﾀｲ ｾｯﾃｲ
20 POKE#700,#08,#23,#1B,#02,#5B,#18,#1A,#68,#04,#33,#18,#68,#10,#60,#58,#68,#50,#60,#98,#68,#90,#60,#D8,#68,#D0,#60,#70,#47
30 LET[0],#0FA8,#1000,0,0,0,0,0,0,0,0
40 FOR I=1 TO 20
50 [8]=I
60 X=USR(#700,0)
70 A=RND(32717)
80 B=RND(32719)
90 C=RND(32749)
100 ?DEC$(I,5);DEC$(A,6);DEC$(B,6);DEC$(C,6)
110 NEXT

```

実行結果は

```
    1     0     0     0
    2     1     1     1
    3     1     1     1
    4     2     2     2
    5     2     2     2
    6     3     3     3
    7     3     3     3
    8     4     4     4
    9     4     4     4
   10     5     5     5
   11     5     5     5
   12     6     6     6
   13     6     6     6
   14     7     7     7
   15     7     7     7
   16     8     8     8
   17     8     8     8
   18     9     9     9
   19     9     9     9
   20    10    10    10
```

中国剰余定理を使うまでもなく、単にXorshiftの結果を2で割った値が「もとの値」になっている…？

40行目を

```
40 FOR I=#1234 TO #1244
```

とし、もう少し大きい値でも実験してみる。結果は

```
 4660  2330  2330  2330
 4661  2330  2330  2330
 4662  2331  2331  2331
 4663  2331  2331  2331
 4664  2332  2332  2332
 4665  2332  2332  2332
 4666  2333  2333  2333
 4667  2333  2333  2333
 4668  2334  2334  2334
 4669  2334  2334  2334
 4670  2335  2335  2335
 4671  2335  2335  2335
 4672  2336  2336  2336
 4673  2336  2336  2336
 4674  2337  2337  2337
 4675  2337  2337  2337
 4676  2338  2338  2338
```

やはり2で割っただけに見える。先ほどの乱数を出力するC言語のプログラムの

```c
for (j = 0; j < 5; j++) printf("%6u", randint() % 0x7fff);
```

という行を、乱数を2で割る

```c
for (j = 0; j < 5; j++) printf("%6u", (randint() / 2) % 0x7fff);
```

に変え、シード`0`で実行すると、結果は

```
  3862 28817 23783 10341 31427
 28768 17208  5212 30730 26758
 10396 28723  6606 28846 10760
 15625 28866  2101  4861 28223
 26797 10169   124  6045 18667
 13713  6134 29440 14295 13830
 22496  2487 19183 13317  8598
 24442 14793 20827  7463 22474
 22358 22888 11557 23760 14514
 12503 25973 19548 17653  8101
 12123 10811 11018 23906 32612
 17319  9595 24906  2513 29018
 10980  4693  9686 27995  9429
 21814  3501 22153 20335 15579
 31676 17659 18705 16419 32000
  1482 16086 11524 28736  4282
 14017   409 18955  6993 20631
 31604 17796 29544 10792 22614
 10946  3863 20139  3208  7284
 24419 19363    51 14622 15347
```

おお！IchigoJamでの実行結果とピッタリ一致した！

# 中途半端な値が出たらやり直す？

ここまでの検証で、IchigoJamでの乱数は「Xorshiftの結果を2で割った値を、`RND(n)`の引数で割った余り」
のようであることがわかった。
しかし、まだ「中途半端な値が出たらやり直す」という処理が入っている可能性がある。
そこで、この可能性を調べる。
おそらくやり直すなら最小値か最大値だろうと思うので、最小値と最大値でやり直すかをチェックする。

最小値の`0`を返させるには、`x`と`w`を`0`に設定しておけばいい。
(Xorshiftの)最大値の`0xffffffff`を返すのは非自明なので、[z3](https://github.com/Z3Prover/z3)を用いて`x`と`w`を求める。

```py
import z3

x = z3.BitVec("x", 32)
w = z3.BitVec("w", 32)

t = x ^ (x << 11)
res = (w ^ z3.LShR(w, 19)) ^ (t ^ z3.LShR(t, 8))

solver = z3.Solver()
solver.add(res == 0xffffffff)

c = solver.check()
print(c)

if c == z3.sat:
	m = solver.model()
	print("x = " + hex(m[x].as_long()))
	print("w = " + hex(m[w].as_long()))
```

結果、`x = 0xc2dc981, w = 0x9dfc4c00`とすればいいようだ。

これらの値を用いてプログラムを組む。
乱数取得後の`0xDEADBEEF`の位置により、Xorshiftの乱数を何個取得したかがわかるはずである。

```
10 ' ｼﾞｮｳﾀｲ ｾｯﾃｲ
20 POKE#700,#08,#23,#1B,#02,#5B,#18,#1A,#68,#04,#33,#18,#68,#10,#60,#58,#68,#50,#60,#98,#68,#90,#60,#D8,#68,#D0,#60,#70,#47
30 ' ｼﾞｮｳﾀｲ ｼｭﾄｸ
40 POKE#780,#08,#23,#1B,#02,#5B,#18,#1A,#68,#04,#33,#10,#68,#18,#60,#50,#68,#58,#60,#90,#68,#98,#60,#D0,#68,#D8,#60,#70,#47
50 LET[0],#0FA8,#1000
60 LET[2],0,0,0,0,#BEEF,#DEAD,0,0
70 X=USR(#700,0)
80 ?RND(100)
90 X=USR(#780,0)
100 ?HEX$([3],4);HEX$([2],4);" ";HEX$([5],4);HEX$([4],4);" ";HEX$([7],4);HEX$([6],4);" ";HEX$([9],4);HEX$([8],4)
110 LET[2],#C981,#0C2D,0,0,#BEEF,#DEAD,#4C00,#9DFC
120 X=USR(#700,0)
130 ?RND(100)
140 X=USR(#780,0)
150 ?HEX$([3],4);HEX$([2],4);" ";HEX$([5],4);HEX$([4],4);" ";HEX$([7],4);HEX$([6],4);" ";HEX$([9],4);HEX$([8],4)
```

実行結果は

```
0
00000000 DEADBEEF 00000000 00000000
47
00000000 DEADBEEF 9DFC4C00 FFFFFFFF
```

となった。きちんと`0`および`0xffffffff`が得られ、Xorshiftは1回ずつしか使われていないことがわかる。
すなわち、このケースではやり直しは起きなかったことがわかる。
`0x7fffffff % 100 == 47`なので、中途半端な値は存在する。
よって、IchigoJamの乱数生成においては、「中途半端な値が出たらやり直す」ことはせず、
1回の乱数生成につき常に1回だけXorshiftを用いていると推測できる。

# 初期化していない時の乱数の状態は？

これで`SRND`で初期化をした後の乱数の挙動はかなり推測できたが、
`SRND`を実行していない時はどうなっているだろうか？
先ほどの乱数の内部状態を取得するプログラムから、`SRND`の行を抜いて実行してみた。

```
20 POKE#700,#08,#23,#1B,#02,#5B,#18,#1A,#68,#04,#33,#10,#68,#18,#60,#50,#68,#58,#60,#90,#68,#98,#60,#D0,#68,#D8,#60,#70,#47
30 LET[0],#0FA8,#1000
50 FOR I=0 TO 20
60 X=USR(#700,0)
70 ?HEX$([3],4);HEX$([2],4);" ";HEX$([5],4);HEX$([4],4);" ";HEX$([7],4);HEX$([6],4);" ";HEX$([9],4);HEX$([8],4)
80 X=RND(100)
90 NEXT
```

実行結果(例)

```
73B23597 9C3E79EE 02AF089C F67EB7CA
9C3E79EE 02AF089C F67EB7CA 14823A1F
02AF089C F67EB7CA 14823A1F 7B1CC068
F67EB7CA 14823A1F 7B1CC068 018DCC7F
14823A1F 7B1CC068 018DCC7F 024EEB63
7B1CC068 018DCC7F 024EEB63 07197BF7
018DCC7F 024EEB63 07197BF7 9A9BE4FC
024EEB63 07197BF7 9A9BE4FC F51A2DE4
07197BF7 9A9BE4FC F51A2DE4 807AD5D7
9A9BE4FC F51A2DE4 807AD5D7 4C70C0EC
F51A2DE4 807AD5D7 4C70C0EC 0989719A
807AD5D7 4C70C0EC 0989719A 2DD80842
4C70C0EC 0989719A 2DD80842 7B5AB443
0989719A 2DD80842 7B5AB443 B1E76C64
2DD80842 7B5AB443 B1E76C64 F3A0DE63
7B5AB443 B1E76C64 F3A0DE63 1ED7424D
B1E76C64 F3A0DE63 1ED7424D B0811578
F3A0DE63 1ED7424D B0811578 3A8FCB40
1ED7424D B0811578 3A8FCB40 CF2959B4
B0811578 3A8FCB40 CF2959B4 6B48AF36
3A8FCB40 CF2959B4 6B48AF36 D3DA5DF2
```

状態の値が左にシフトしていき、初期化した時と同様にXorshiftの計算が1回ずつ行われていると推測できる。
まあ、ここをわざわざ変えるメリットは薄いだろう。
さらにループも外して、実行開始時の乱数の状態をいくつか見てみる。

```
10 POKE#700,#08,#23,#1B,#02,#5B,#18,#1A,#68,#04,#33,#10,#68,#18,#60,#50,#68,#58,#60,#90,#68,#98,#60,#D0,#68,#D8,#60,#70,#47
20 LET[0],#0FA8,#1000
30 X=USR(#700,0)
40 ?HEX$([3],4);HEX$([2],4);" ";HEX$([5],4);HEX$([4],4);" ";HEX$([7],4);HEX$([6],4);" ";HEX$([9],4);HEX$([8],4)
```

実行結果(例)

```
RUN
D64A814C EEDCE49A 03741E59 6C556546
OK
RUN
47897E28 E77C5423 BFBC360D 0A6E336F
OK
RUN
531E5DCB 66F54604 8A92A23C 62378948
OK
RUN
2572921B E69CD6A3 0DEF975E 3F21DE2F
OK
RUN
A865C84A D9B762D6 B3D4D9D3 39910975
OK
RUN
805E61A5 10BCF767 C4BCD6DD CD16535D
OK
RUN
506914C7 4E31C623 E68B4B26 62674095
OK
```

うーん…バラバラの値になっており、法則はよくわからない。

# 結論

今回は、以下のことが推測できた。

* IchigoJamの`RND(n)`関数は、「Xorshiftの乱数を2で割った値を、`n`で割った余り」を返す。
  * `n`が0以下の時は、かわりに0を返す。
  * `SRND n`命令が無い古いバージョンについては不明。
* Xorshiftから極端に大きい値や極端に小さい値が返されても、やり直さずそのまま使う。
* `SRND n`は、Xorshiftの状態を`x = n, y = 362436069, z = 521288629, w = 88675123`に設定する。
  * `x`は符号拡張される。
* インタラクティブからプログラムを実行開始する時、Xorshiftの状態が適当な値に設定される。
  * 特に特定の部分を固定する様子はみられない。

以前[IchigoJamのマシン語とBASICでXorshiftを実装したことがあった](https://qiita.com/mikecat_mixc/items/e674b6d21f9df4c9f730)が、
IchigoJam自体でもXorshiftが使われているっぽいとはなあ…。
それだけ便利なアルゴリズムなんだなあ…。

※IchigoJamはjig.jpの登録商標です。
※IchigoKamuyは株式会社syushuの登録商標です。