[閒聊] 替RandomVariate提供提取seed的功能

作者: AmibaGelos (Amiba Gelos)   2016-10-11 12:28:18
最近在用MMA跑Monte Carlo,覺得MMA不能提取seed實在是很不方便,所以就做了一個可以
提取seed的RandomVariate
Module[{seed, store = {}, temp = 1, Storing},
RandomMemory[1, opt___] := RandomMemory[temp, opt];
RandomMemory[in__] := (Storing[in]; RandomVariate[in]);
Storing[dist_, opt___] :=
AppendTo[
store, {Evaluate[If[dist === temp, 1, temp = dist]], opt}];
SetSeed[InSeed_, InStore_: {}] := (store = {}; temp = 1;
SeedRandom[seed = InSeed]; (RandomMemory[##];) & @@@ InStore;);
ExtractSeed[] := {seed, store}; SetSeed[0];];
這個設計在跑custom distribution的時候蠻好用的, 不過AppendTo會是瓶頸
目前只有想到hash, 不過有新distribution還是得重建table
不知道有沒有更漂亮的寫法?
作者: ginstein (邁向學術之路)   2016-10-11 21:30:00
不懂! seed不就使用者設定的嗎? 不能提取自己設的seed?不懂!seed不就使用者設定的嗎?不能提取自己設的seed?
作者: LPH66 (-6.2598534e+18f)   2016-10-12 03:29:00
他要的是跑到一半記錄目前 PRNG 的 state 吧不過記得不少演算法的 state 並不只有 seed 所以有點難提供一個統一界面, 例如 MT19937 需要 624 個 32 位元整數但種子只有 32 位元, 所以不可能用單一整數表示狀態然後, 剛剛翻 Help 有看到 BlockRandom 或許也符合你的需要
作者: AmibaGelos (Amiba Gelos)   2016-10-12 09:57:00
謝謝樓上解釋~MMA預設CA是80*64b,最小是64*64所以不可能用long描述, MMA的seed我覺得實際上是被hash處理過這方面我還在測試,不過想從這邊retrive state也是很難BlockRandom只是freeze state而已,應該不能做到同時跑兩個prng,更做不到resume吧@@其實有種最偷懶的做法是每次要存檔的時候就另開一個kernel專門跑rng,用mathlink連到另一個MMA也許更clear忽然覺得這應該就是最佳解了@@唯一可惜的是還是得額外存一套distribution,跑MCMC有點浪費ram...
作者: ginstein (邁向學術之路)   2016-10-12 18:31:00
謝謝L兄,懂了.另MMAseed代表所有狀態,沒看過找seed指令

Links booklink

Contact Us: admin [ a t ] ucptt.com