auto merge of #10196 : huonw/rust/fix-zig, r=alexcrichton

The code was using (in the notation of Doornik 2005) `f(x_{i+1}) -
f(x_{i+2})` rather than `f(x_i) - f(x_{i+1})`. This corrects that, and
removes the F_DIFF tables which caused this problem in the first place.

They `F_DIFF` tables are a micro-optimisation (in theory, they could
easily be a micro-pessimisation): that `if` gets hit about 1% of the
time for Exp/Normal, and the rest of the condition involves RNG calls
and a floating point `exp`, so it is unlikely that saving a single FP
subtraction will be very useful (especially as more tables means more
memory reads and higher cache pressure, as well as taking up space in
the binary (although only ~2k in this case)).

Closes #10084. Notably, unlike that issue suggests, this wasn't a
problem with the Exp tables. It affected Normal too, but since it is
symmetric, there was no bias in the mean (as the bias was equal on the
positive and negative sides and so cancelled out) but it was visible as
a variance slightly lower than it should be.

New plot:

![exp-density](https://f.cloud.github.com/assets/1203825/1445796/42218dfe-422a-11e3-9f98-2cd146b82b46.png)

I've started writing some tests in [huonw/random-tests](https://github.com/huonw/random-tests) (not in the main repo because they can and do fail occasionally, due to randomness, but it is on Travis and Rust-CI so it will hopefully track the language), unsurprisingly, they're [currently failing](https://travis-ci.org/huonw/random-tests/builds/13313987) (note that both exp and norm are failing, the former due to both mean and variance the latter due to just variance), but pass at the 0.01 level reliably with this change.

(Currently the only test is essentially a quantitative version of the plots I've been showing, which is run on the `f64` `Rand` instance (uniform 0 to 1), and the Normal and Exp distributions.)

diff --git a/src/etc/ziggurat_tables.py b/src/etc/ziggurat_tables.py
index 51c0da39bd5ac167df3222ce41c469389870b137..263c76a88ddcab885664cbe9452dcf55a6293003 100755 (executable)
--- a/src/etc/ziggurat_tables.py
+++ b/src/etc/ziggurat_tables.py
@@ -19,7 +19,7 @@ from math import exp, sqrt, log, floor
 import random
 
 # The order should match the return value of `tables`
-TABLE_NAMES = ['X', 'F', 'F_DIFF']
+TABLE_NAMES = ['X', 'F']
 
 # The actual length of the table is 1 more, to stop
 # index-out-of-bounds errors. This should match the bitwise operation
@@ -43,13 +43,10 @@ def tables(r, v, f, f_inv):
 
     # cache the f's
     fvec = [0]*(TABLE_LEN+1)
-    fdiff = [0]*(TABLE_LEN+1)
     for i in range(TABLE_LEN+1):
         fvec[i] = f(xvec[i])
-        if i > 0:
-            fdiff[i] = fvec[i] - fvec[i-1]
 
-    return xvec, fvec, fdiff
+    return xvec, fvec
 
 # Distributions
 # N(0, 1)

diff --git a/src/libstd/rand/distributions.rs b/src/libstd/rand/distributions.rs
index e7bcf8ce5d3a71d0c0d14e41b6199cc0a1093c3a..a6fa1a2bd875fb3220d148ba36aedd15d75103b7 100644 (file)
--- a/src/libstd/rand/distributions.rs
+++ b/src/libstd/rand/distributions.rs
@@ -207,7 +207,6 @@ fn ziggurat<R:Rng>(rng: &mut R,
                    symmetric: bool,
                    X: ziggurat_tables::ZigTable,
                    F: ziggurat_tables::ZigTable,
-                   F_DIFF: ziggurat_tables::ZigTable,
                    pdf: &'static fn(f64) -> f64,
                    zero_case: &'static fn(&mut R, f64) -> f64) -> f64 {
     static SCALE: f64 = (1u64 << 53) as f64;
@@ -237,7 +236,7 @@ fn ziggurat<R:Rng>(rng: &mut R,
             return zero_case(rng, u);
         }
         // algebraically equivalent to f1 + DRanU()*(f0 - f1) < 1
-        if F[i + 1] + F_DIFF[i + 1] * rng.gen() < pdf(x) {
+        if F[i + 1] + (F[i] - F[i + 1]) * rng.gen() < pdf(x) {
             return x;
         }
     }
@@ -288,7 +287,7 @@ fn zero_case<R:Rng>(rng: &mut R, u: f64) -> f64 {
             rng,
             true, // this is symmetric
             &ziggurat_tables::ZIG_NORM_X,
-            &ziggurat_tables::ZIG_NORM_F, &ziggurat_tables::ZIG_NORM_F_DIFF,
+            &ziggurat_tables::ZIG_NORM_F,
             pdf, zero_case))
     }
 }
@@ -366,7 +365,7 @@ fn zero_case<R:Rng>(rng: &mut R, _u: f64) -> f64 {
 
         Exp1(ziggurat(rng, false,
                       &ziggurat_tables::ZIG_EXP_X,
-                      &ziggurat_tables::ZIG_EXP_F, &ziggurat_tables::ZIG_EXP_F_DIFF,
+                      &ziggurat_tables::ZIG_EXP_F,
                       pdf, zero_case))
     }
 }

diff --git a/src/libstd/rand/ziggurat_tables.rs b/src/libstd/rand/ziggurat_tables.rs
index aca2457cac42c282dc10355e1ba38e827a7da40f..049ef3dbb5936933b894b9dfe9891fa872affa97 100644 (file)
--- a/src/libstd/rand/ziggurat_tables.rs
+++ b/src/libstd/rand/ziggurat_tables.rs
@@ -145,72 +145,6 @@
      0.887984660763399880, 0.898095921906304051, 0.908726440060562912, 0.919991505048360247,
      0.932060075968990209, 0.945198953453078028, 0.959879091812415930, 0.977101701282731328,
      1.000000000000000000];
-pub static ZIG_NORM_F_DIFF: [f64, .. 257] =
-    [0.000000000000000000, 0.000782818165911943, 0.001348786815607765, 0.001428899847265509,
-     0.001484430705892882, 0.001528472172127356, 0.001565707298030807, 0.001598388670308183,
-     0.001627786418212004, 0.001654692743837703, 0.001679637706201265, 0.001702994844613767,
-     0.001725038123187510, 0.001745974954326004, 0.001765966477270568, 0.001785140598493315,
-     0.001803600702759419, 0.001821431661060659, 0.001838704088536796, 0.001855477433793579,
-     0.001871802266665008, 0.001887722003144375, 0.001903274226858077, 0.001918491715965767,
-     0.001933403251421835, 0.001948034260540625, 0.001962407334827158, 0.001976542650643127,
-     0.001990458313945481, 0.002004170645086643, 0.002017694415851860, 0.002031043048104267,
-     0.002044228781321551, 0.002057262814738517, 0.002070155428613822, 0.002082916088226049,
-     0.002095553533492583, 0.002108075856553551, 0.002120490569226280, 0.002132804661891696,
-     0.002145024655099026, 0.002157156644953973, 0.002169206343177243, 0.002181179112575302,
-     0.002193079998548175, 0.002204913757158977, 0.002216684880213121, 0.002228397617726446,
-     0.002240055998106505, 0.002251663846325885, 0.002263224800326716, 0.002274742325862292,
-     0.002286219729956393, 0.002297660173134250, 0.002309066680560787, 0.002320442152205823,
-     0.002331789372137141, 0.002343111017035562, 0.002354409664009627, 0.002365687797781804,
-     0.002376947817308683, 0.002388192041889739, 0.002399422716815966, 0.002410642018598946,
-     0.002421852059823287, 0.002433054893654529, 0.002444252518034679, 0.002455446879594508,
-     0.002466639877306970, 0.002477833365903986, 0.002489029159078809, 0.002500229032490808,
-     0.002511434726590794, 0.002522647949281448, 0.002533870378427505, 0.002545103664226889,
-     0.002556349431455662, 0.002567609281597438, 0.002578884794865288, 0.002590177532127119,
-     0.002601489036740262, 0.002612820836305291, 0.002624174444343735, 0.002635551361907296,
-     0.002646953079123743, 0.002658381076686089, 0.002669836827288052, 0.002681321797012387,
-     0.002692837446676144, 0.002704385233135737, 0.002715966610556786, 0.002727583031652520,
-     0.002739235948893221, 0.002750926815690169, 0.002762657087557796, 0.002774428223256353,
-     0.002786241685917290, 0.002798098944155558, 0.002810001473169871, 0.002821950755833219,
-     0.002833948283778004, 0.002845995558475284, 0.002858094092312607, 0.002870245409671041,
-     0.002882451048004164, 0.002894712558920987, 0.002907031509275432, 0.002919409482262880,
-     0.002931848078526783, 0.002944348917277934, 0.002956913637427061, 0.002969543898733384,
-     0.002982241382970874, 0.002995007795115689, 0.003007844864553855, 0.003020754346314269,
-     0.003033738022328147, 0.003046797702715820, 0.003059935227105459, 0.003073152465984053,
-     0.003086451322084072, 0.003099833731808721, 0.003113301666695822, 0.003126857134927052,
-     0.003140502182881588, 0.003154238896738770, 0.003168069404132778, 0.003181995875862154,
-     0.003196020527657495, 0.003210145622009941, 0.003224373470066433, 0.003238706433592253,
-     0.003253146927007733, 0.003267697419501892, 0.003282360437226572, 0.003297138565578506,
-     0.003312034451571411, 0.003327050806304299, 0.003342190407532641, 0.003357456102345890,
-     0.003372850809960137, 0.003388377524629727, 0.003404039318688046, 0.003419839345721265,
-     0.003435780843885239, 0.003451867139373843, 0.003468101650046629, 0.003484487889225119,
-     0.003501029469670069, 0.003517730107746697, 0.003534593627793237, 0.003551623966702611,
-     0.003568825178730639, 0.003586201440546166, 0.003603757056536316, 0.003621496464384588,
-     0.003639424240937217, 0.003657545108379068, 0.003675863940735269, 0.003694385770723563,
-     0.003713115796977806, 0.003732059391668707, 0.003751222108547281, 0.003770609691440940,
-     0.003790228083232539, 0.003810083435355216, 0.003830182117840641, 0.003850530729957835,
-     0.003871136111486317, 0.003892005354668437, 0.003913145816891062, 0.003934565134149914,
-     0.003956271235355358, 0.003978272357543333, 0.004000577062061084, 0.004023194251800533,
-     0.004046133189565926, 0.004069403517661885, 0.004093015278800460, 0.004116978938436600,
-     0.004141305408647655, 0.004166006073685835, 0.004191092817346642, 0.004216578052307351,
-     0.004242474751606884, 0.004268796482457593, 0.004295557442594244, 0.004322772499391836,
-     0.004350457232007221, 0.004378627976825644, 0.004407301876525049, 0.004436496933105327,
-     0.004466232065271192, 0.004496527170598785, 0.004527403192966406, 0.004558882195791591,
-     0.004590987441673855, 0.004623743479123199, 0.004657176237135574, 0.004691313128472929,
-     0.004726183162616859, 0.004761817069491636, 0.004798247435199299, 0.004835508851176451,
-     0.004873638078381815, 0.004912674228345848, 0.004952658963181422, 0.004993636716962402,
-     0.005035654941235035, 0.005078764377854039, 0.005123019362831771, 0.005168478165478940,
-     0.005215203367812893, 0.005263262290042703, 0.005312727468930079, 0.005363677197016692,
-     0.005416196132139284, 0.005470375988385734, 0.005526316321746716, 0.005584125426278286,
-     0.005643921359735682, 0.005705833121505521, 0.005770002010457520, 0.005836583196307310,
-     0.005905747545561058, 0.005977683752542928, 0.006052600837980204, 0.006130731092920838,
-     0.006212333565464245, 0.006297698213369562, 0.006387150879090475, 0.006481059288027780,
-     0.006579840329791975, 0.006683968961788356, 0.006793989182803495, 0.006910527673723577,
-     0.007034310911336661, 0.007166186857056056, 0.007307152748134871, 0.007458391141830445,
-     0.007621317291194862, 0.007797642342679434, 0.007989459040836144, 0.008199360125510702,
-     0.008430605346682607, 0.008687362737884952, 0.008975066840784529, 0.009300967772353674,
-     0.009675004947253041, 0.010111261142904171, 0.010630518154258861, 0.011265064987797335,
-     0.012068570920629962, 0.013138877484087819, 0.014680138359337902, 0.017222609470315398,
-     0.022898298717268672];
 pub static ZIG_EXP_R: f64 = 7.697117470131050077;
 pub static ZIG_EXP_X: [f64, .. 257] =
     [8.697117470131052741, 7.697117470131050077, 6.941033629377212577, 6.478378493832569696,
@@ -344,69 +278,3 @@
      0.775956852040116218, 0.791527636972496285, 0.808421651523009044, 0.826993296643051101,
      0.847785500623990496, 0.871704332381204705, 0.900469929925747703, 0.938143680862176477,
      1.000000000000000000];
-pub static ZIG_EXP_F_DIFF: [f64, .. 257] =
-    [0.000000000000000000, 0.000287067661533533, 0.000513134928485678, 0.000569030497974398,
-     0.000609667963417335, 0.000642831049855169, 0.000671465984262828, 0.000697030342996893,
-     0.000720360862708599, 0.000741986223663093, 0.000762263730113694, 0.000781447246315807,
-     0.000799724254382053, 0.000817237547791934, 0.000834098656693235, 0.000850396538527769,
-     0.000866203416804620, 0.000881578828420777, 0.000896572504999613, 0.000911226471926952,
-     0.000925576608509206, 0.000939653828282008, 0.000953484986066785, 0.000967093584871414,
-     0.000980500333784669, 0.000993723593313716, 0.001006779734568374, 0.001019683431705467,
-     0.001032447902101660, 0.001045085105172934, 0.001057605908173612, 0.001070020225402434,
-     0.001082337135821582, 0.001094564983022843, 0.001106711460658764, 0.001118783685829211,
-     0.001130788262427001, 0.001142731336065933, 0.001154618641914802, 0.001166455546523074,
-     0.001178247084534012, 0.001189997991027938, 0.001201712730115490, 0.001213395520299268,
-     0.001225050357040701, 0.001236681032901414, 0.001248291155571943, 0.001259884164055092,
-     0.001271463343231895, 0.001283031837006378, 0.001294592660197942, 0.001306148709326875,
-     0.001317702772419903, 0.001329257537945404, 0.001340815602974395, 0.001352379480650950,
-     0.001363951607045839, 0.001375534347457789, 0.001387130002219621, 0.001398740812059381,
-     0.001410368963061376, 0.001422016591266340, 0.001433685786946429, 0.001445378598586011,
-     0.001457097036596827, 0.001468843076792140, 0.001480618663643060, 0.001492425713336909,
-     0.001504266116655995, 0.001516141741693663, 0.001528054436422108, 0.001540006031125918,
-     0.001551998340713470, 0.001564033166917514, 0.001576112300394977, 0.001588237522735750,
-     0.001600410608388780, 0.001612633326513305, 0.001624907442762655, 0.001637234721007311,
-     0.001649616925003372, 0.001662055820012304, 0.001674553174376953, 0.001687110761059388,
-     0.001699730359144919, 0.001712413755316500, 0.001725162745304071, 0.001737979135312442,
-     0.001750864743431488, 0.001763821401032123, 0.001776850954151601, 0.001789955264870927,
-     0.001803136212688003, 0.001816395695889220, 0.001829735632922019, 0.001843157963772116,
-     0.001856664651347151, 0.001870257682870316, 0.001883939071285826, 0.001897710856679738,
-     0.001911575107717528, 0.001925533923102574, 0.001939589433056721, 0.001953743800826108,
-     0.001967999224215228, 0.001982357937151347, 0.001996822211282223, 0.002011394357609747,
-     0.002026076728162574, 0.002040871717710169, 0.002055781765521847, 0.002070809357173103,
-     0.002085957026402963, 0.002101227357025226, 0.002116622984897121, 0.002132146599948981,
-     0.002147800948277823, 0.002163588834309782, 0.002179513123034188, 0.002195576742314159,
-     0.002211782685277469, 0.002228134012792427, 0.002244633856033434, 0.002261285419141418,
-     0.002278091981983449, 0.002295056903017983, 0.002312183622271174, 0.002329475664429648,
-     0.002346936642057179, 0.002364570258941101, 0.002382380313575932, 0.002400370702791893,
-     0.002418545425535629, 0.002436908586812392, 0.002455464401797752, 0.002474217200128692,
-     0.002493171430384328, 0.002512331664766249, 0.002531702603989994, 0.002551289082400404,
-     0.002571096073321844, 0.002591128694658967, 0.002611392214760672, 0.002631892058563845,
-     0.002652633814032662, 0.002673623238910738, 0.002694866267805934, 0.002716369019626269,
-     0.002738137805389534, 0.002760179136428037, 0.002782499733014893, 0.002805106533435520,
-     0.002828006703534697, 0.002851207646767162, 0.002874717014785921, 0.002898542718600849,
-     0.002922692940346749, 0.002947176145699226, 0.002972001096982591, 0.002997176867015228,
-     0.003022712853742948, 0.003048618795714386, 0.003074904788455568, 0.003101581301807876,
-     0.003128659198296080, 0.003156149752600867, 0.003184064672214937, 0.003212416119368622,
-     0.003241216734320596, 0.003270479660111680, 0.003300218568896729, 0.003330447689969929,
-     0.003361181839619420, 0.003392436452949343, 0.003424227617828290, 0.003456572111131984,
-     0.003489487437467131, 0.003522991870580083, 0.003557104497672658, 0.003591845266868621,
-     0.003627235038102472, 0.003663295637722386, 0.003700049917134574, 0.003737521815846301,
-     0.003775736429304177, 0.003814720081962375, 0.003854500406067995, 0.003895106426696382,
-     0.003936568653631844, 0.003978919180756157, 0.004022191793678687, 0.004066422086428989,
-     0.004111647588127876, 0.004157907900659452, 0.004205244848493050, 0.004253702641940915,
-     0.004303328055299205, 0.004354170621502118, 0.004406282845128784, 0.004459720435841752,
-     0.004514542564613699, 0.004570812145417769, 0.004628596145424491, 0.004687965927177740,
-     0.004748997626717266, 0.004811772572194672, 0.004876377748206484, 0.004942906311860507,
-     0.005011458167522187, 0.005082140608288488, 0.005155069033533799, 0.005230367753417398,
-     0.005308170893076836, 0.005388623411430704, 0.005471882252147620, 0.005558117647517014,
-     0.005647514599798176, 0.005740274569295156, 0.005836617404105682, 0.005936783553485037,
-     0.006041036615386131, 0.006149666279423593, 0.006262991739818591, 0.006381365669577810,
-     0.006505178868201678, 0.006634865721946159, 0.006770910649812723, 0.006913855752425535,
-     0.007064309938019209, 0.007222959874423007, 0.007390583214465396, 0.007568064673498798,
-     0.007756415714389786, 0.007956798835585532, 0.008170557788458321, 0.008399255510700199,
-     0.008644722212900025, 0.008909116987305010, 0.009195007664428712, 0.009505475652925033,
-     0.009844255532840629, 0.010215923852312625, 0.010626158965710175, 0.011082105722287849,
-     0.011592898788496009, 0.012170432837851575, 0.012830529553771619, 0.013594766864701180,
-     0.014493463190219380, 0.015570784932380066, 0.016894014550512759, 0.018571645120042057,
-     0.020792203980939394, 0.023918831757214210, 0.028765597544542998, 0.037673750936428774,
-     0.061856319137823523];