Last active
September 29, 2025 13:43
-
-
Save rnapier/f70ad0beb0780f974fbac7267e297508 to your computer and use it in GitHub Desktop.
Ziggurat algorithm ported from Java using Codex: https://mastodon.social/@cocoaphony/115282998859062341
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| /// Test to check sample distribution. Written by Codex. | |
| /// Prompt: | |
| /// Write unit tests for NormalizedGaussian with the goal of checking that it generates random values with a guassian normal distribution. | |
| import XCTest | |
| @testable import FinLib | |
| final class NormalizedGaussianTests: XCTestCase { | |
| private let sampleCount = 200_000 | |
| func testMeanAndVarianceMatchStandardNormal() { | |
| var sampler = makeSampler(seed: 0x1234_5678_9ABC_DEF0) | |
| var sum = 0.0 | |
| var sumSquares = 0.0 | |
| for _ in 0..<sampleCount { | |
| let value = sampler.sample() | |
| sum += value | |
| sumSquares += value * value | |
| } | |
| let mean = sum / Double(sampleCount) | |
| let variance = sumSquares / Double(sampleCount) - mean * mean | |
| XCTAssertEqual(mean, 0, accuracy: 0.01) | |
| XCTAssertEqual(variance, 1, accuracy: 0.02) | |
| } | |
| func testEmpiricalRuleCoverage() { | |
| var sampler = makeSampler(seed: 0x0FED_CBA9_8765_4321) | |
| var withinOne = 0 | |
| var withinTwo = 0 | |
| var withinThree = 0 | |
| for _ in 0..<sampleCount { | |
| let value = sampler.sample() | |
| let magnitude = abs(value) | |
| if magnitude <= 1 { withinOne += 1 } | |
| if magnitude <= 2 { withinTwo += 1 } | |
| if magnitude <= 3 { withinThree += 1 } | |
| } | |
| let total = Double(sampleCount) | |
| XCTAssertEqual(Double(withinOne) / total, 0.6827, accuracy: 0.01) | |
| XCTAssertEqual(Double(withinTwo) / total, 0.9545, accuracy: 0.01) | |
| XCTAssertEqual(Double(withinThree) / total, 0.9973, accuracy: 0.005) | |
| } | |
| private func makeSampler(seed: UInt64) -> ZigguratSamplers.NormalizedGaussian { | |
| ZigguratSamplers.NormalizedGaussian(rng: SplitMix64Random(seed: seed)) | |
| } | |
| } | |
| private struct SplitMix64Random: UniformRandomProvider { | |
| private var state: UInt64 | |
| init(seed: UInt64) { | |
| state = seed | |
| } | |
| mutating func nextLong() -> UInt64 { | |
| state = state &+ 0x9E3779B97F4A7C15 | |
| var z = state | |
| z = (z ^ (z >> 30)) &* 0xBF58476D1CE4E5B9 | |
| z = (z ^ (z >> 27)) &* 0x94D049BB133111EB | |
| return z ^ (z >> 31) | |
| } | |
| } |
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| /// Based closely on https://commons.apache.org/proper/commons-rng/commons-rng-docs/jacoco-aggregate/commons-rng-sampling/org.apache.commons.rng.sampling.distribution/ZigguratSampler.java.html | |
| /// For details on Ziggurat: https://en.wikipedia.org/wiki/Ziggurat_algorithm | |
| /// This provided about a 7% performance boost over my previous implementation of Marsaglia Polar (which was a big improvement over Box Muller). | |
| /// | |
| /// Written by Codex. | |
| /// Primary Prompt: | |
| /// Convert ZigguratSampler.java to Swift. Keep the implementation as close to the original Java as possible. If you find anything missing, ask me about it. | |
| /// | |
| /// Later changes using Codex included: | |
| /// * converting classes to structs | |
| /// * removing public annotations (to get better inlining) | |
| /// * removing unused methods (done by hand, not Codex) | |
| /// * Converted from Int64 to UInt64 handling. (Getting this right by hand was tricky; Codex did a much better job at getting it right.) | |
| /// * Rename nextLong() to nextUniform() to match other code | |
| import Foundation | |
| /// Protocol defining a source of 64-bit uniformly distributed random numbers. | |
| protocol UniformRandomProvider { | |
| mutating func nextUniform() -> UInt64 | |
| } | |
| /// Shared interface for samplers that rely on the ziggurat algorithm. | |
| protocol ZigguratSampler { | |
| var rng: UniformRandomProvider { get set } | |
| } | |
| extension ZigguratSampler { | |
| /// Mask to extract the lowest 8-bits from an integer. | |
| static var maskInt8: Int { 0xff } | |
| /// 2^63. | |
| static var twoPow63: Double { Double(0x1.0p63) } | |
| mutating func randomInt63() -> UInt64 { | |
| rng.nextUniform() >> 1 | |
| } | |
| static func interpolate(_ v: [Double], j: Int, u: UInt64) -> Double { | |
| return v[j] * twoPow63 + Double(u) * (v[j - 1] - v[j]) | |
| } | |
| } | |
| /// Namespace for concrete ziggurat-based samplers. | |
| enum ZigguratSamplers {} | |
| extension ZigguratSamplers { | |
| /// Modified ziggurat method for sampling from an exponential distribution. | |
| struct Exponential: ZigguratSampler { | |
| /// Number of layers in the ziggurat. Maximum i value for early exit. | |
| private static let I_MAX = 252 | |
| /// Constant used when sampling the tail. | |
| private static let X_0 = 7.569274694148063 | |
| /// Maximum epsilon distance of pdf(x) from the triangle hypotenuse for early acceptance. | |
| private static let E_MAX: UInt64 = 853965788476313647 | |
| /// The precomputed ziggurat lengths, denoted X_i in the main text. | |
| private static let X: [Double] = [ | |
| /* [ 0] */ 8.2066240675348816e-19, 7.3973732351607284e-19, 6.9133313377915293e-19, 6.5647358820964533e-19, | |
| /* [ 4] */ 6.2912539959818508e-19, 6.0657224129604964e-19, 5.8735276103737269e-19, 5.7058850528536941e-19, | |
| /* [ 8] */ 5.557094569162239e-19, 5.4232438903743953e-19, 5.3015297696508776e-19, 5.1898739257708062e-19, | |
| /* [ 12] */ 5.086692261799833e-19, 4.9907492938796469e-19, 4.9010625894449536e-19, 4.8168379010649187e-19, | |
| /* [ 16] */ 4.7374238653644714e-19, 4.6622795807196824e-19, 4.5909509017784048e-19, 4.5230527790658154e-19, | |
| /* [ 20] */ 4.458255881635396e-19, 4.3962763126368381e-19, 4.336867596710647e-19, 4.2798143618469714e-19, | |
| /* [ 24] */ 4.2249273027064889e-19, 4.172039125346411e-19, 4.1210012522465616e-19, 4.0716811225869233e-19, | |
| /* [ 28] */ 4.0239599631006903e-19, 3.9777309342877357e-19, 3.9328975785334499e-19, 3.8893725129310323e-19, | |
| /* [ 32] */ 3.8470763218720385e-19, 3.8059366138180143e-19, 3.765887213854473e-19, 3.7268674692030177e-19, | |
| /* [ 36] */ 3.6888216492248162e-19, 3.6516984248800068e-19, 3.6154504153287473e-19, 3.5800337915318032e-19, | |
| /* [ 40] */ 3.5454079284533432e-19, 3.5115350988784242e-19, 3.4783802030030962e-19, 3.4459105288907336e-19, | |
| /* [ 44] */ 3.4140955396563316e-19, 3.3829066838741162e-19, 3.3523172262289001e-19, 3.3223020958685874e-19, | |
| /* [ 48] */ 3.2928377502804472e-19, 3.2639020528202049e-19, 3.2354741622810815e-19, 3.2075344331080789e-19, | |
| /* [ 52] */ 3.1800643250478609e-19, 3.1530463211820845e-19, 3.1264638534265134e-19, 3.1003012346934211e-19, | |
| /* [ 56] */ 3.0745435970137301e-19, 3.0491768350005559e-19, 3.0241875541094565e-19, 2.999563023214455e-19, | |
| /* [ 60] */ 2.9752911310742592e-19, 2.9513603463113224e-19, 2.9277596805684267e-19, 2.9044786545442563e-19, | |
| /* [ 64] */ 2.8815072666416712e-19, 2.8588359639906928e-19, 2.8364556156331615e-19, 2.8143574876779799e-19, | |
| /* [ 68] */ 2.7925332202553125e-19, 2.7709748061152879e-19, 2.7496745707320232e-19, 2.7286251537873397e-19, | |
| /* [ 72] */ 2.7078194919206054e-19, 2.687250802641905e-19, 2.6669125693153442e-19, 2.6467985271278891e-19, | |
| /* [ 76] */ 2.6269026499668434e-19, 2.6072191381359757e-19, 2.5877424068465143e-19, 2.5684670754248168e-19, | |
| /* [ 80] */ 2.5493879571835479e-19, 2.5305000499077481e-19, 2.511798526911271e-19, 2.4932787286227806e-19, | |
| /* [ 84] */ 2.474936154663866e-19, 2.4567664563848669e-19, 2.4387654298267842e-19, 2.4209290090801527e-19, | |
| /* [ 88] */ 2.4032532600140538e-19, 2.3857343743505147e-19, 2.3683686640614648e-19, 2.3511525560671253e-19, | |
| /* [ 92] */ 2.3340825872163284e-19, 2.3171553995306794e-19, 2.3003677356958333e-19, 2.2837164347843482e-19, | |
| /* [ 96] */ 2.2671984281957174e-19, 2.2508107358001938e-19, 2.2345504622739592e-19, 2.2184147936140775e-19, | |
| /* [100] */ 2.2024009938224424e-19, 2.1865064017486842e-19, 2.1707284280826716e-19, 2.1550645524878675e-19, | |
| /* [104] */ 2.1395123208673778e-19, 2.124069342755064e-19, 2.1087332888245875e-19, 2.0935018885097035e-19, | |
| /* [108] */ 2.0783729277295508e-19, 2.0633442467130712e-19, 2.0484137379170616e-19, 2.0335793440326865e-19, | |
| /* [112] */ 2.018839056075609e-19, 2.0041909115551697e-19, 1.9896329927183254e-19, 1.975163424864309e-19, | |
| /* [116] */ 1.9607803747261946e-19, 1.9464820489157862e-19, 1.9322666924284314e-19, 1.9181325872045647e-19, | |
| /* [120] */ 1.9040780507449479e-19, 1.8901014347767504e-19, 1.8762011239677479e-19, 1.8623755346860768e-19, | |
| /* [124] */ 1.8486231138030984e-19, 1.8349423375370566e-19, 1.8213317103353295e-19, 1.8077897637931708e-19, | |
| /* [128] */ 1.7943150556069476e-19, 1.7809061685599652e-19, 1.7675617095390567e-19, 1.7542803085801941e-19, | |
| /* [132] */ 1.7410606179414531e-19, 1.727901311201724e-19, 1.7148010823836362e-19, 1.7017586450992059e-19, | |
| /* [136] */ 1.6887727317167824e-19, 1.6758420925479093e-19, 1.6629654950527621e-19, 1.6501417230628659e-19, | |
| /* [140] */ 1.6373695760198277e-19, 1.624647868228856e-19, 1.6119754281258616e-19, 1.5993510975569615e-19, | |
| /* [144] */ 1.5867737310692309e-19, 1.5742421952115544e-19, 1.5617553678444595e-19, 1.5493121374578016e-19, | |
| /* [148] */ 1.5369114024951992e-19, 1.5245520706841019e-19, 1.5122330583703858e-19, 1.4999532898563561e-19, | |
| /* [152] */ 1.4877116967410352e-19, 1.4755072172615974e-19, 1.4633387956347966e-19, 1.4512053813972103e-19, | |
| /* [156] */ 1.4391059287430991e-19, 1.4270393958586506e-19, 1.4150047442513381e-19, 1.4030009380730888e-19, | |
| /* [160] */ 1.3910269434359025e-19, 1.3790817277185197e-19, 1.3671642588626657e-19, 1.3552735046573446e-19, | |
| /* [164] */ 1.3434084320095729e-19, 1.3315680061998685e-19, 1.3197511901207148e-19, 1.3079569434961214e-19, | |
| /* [168] */ 1.2961842220802957e-19, 1.2844319768333099e-19, 1.2726991530715219e-19, 1.2609846895903523e-19, | |
| /* [172] */ 1.2492875177568625e-19, 1.237606560569394e-19, 1.2259407316813331e-19, 1.2142889343858445e-19, | |
| /* [176] */ 1.2026500605581765e-19, 1.1910229895518744e-19, 1.1794065870449425e-19, 1.1677997038316715e-19, | |
| /* [180] */ 1.1562011745554883e-19, 1.1446098163777869e-19, 1.1330244275772562e-19, 1.1214437860737343e-19, | |
| /* [184] */ 1.109866647870073e-19, 1.0982917454048923e-19, 1.0867177858084351e-19, 1.0751434490529747e-19, | |
| /* [188] */ 1.0635673859884002e-19, 1.0519882162526621e-19, 1.0404045260457141e-19, 1.0288148657544097e-19, | |
| /* [192] */ 1.0172177474144965e-19, 1.0056116419943559e-19, 9.9399497648346677e-20, 9.8236613076667446e-20, | |
| /* [196] */ 9.7072343426320094e-20, 9.5906516230690634e-20, 9.4738953224154196e-20, 9.3569469920159036e-20, | |
| /* [200] */ 9.2397875154569468e-20, 9.1223970590556472e-20, 9.0047550180852874e-20, 8.8868399582647627e-20, | |
| /* [204] */ 8.768629551976745e-20, 8.6501005086071005e-20, 8.5312284983141187e-20, 8.4119880684385214e-20, | |
| /* [208] */ 8.292352551651342e-20, 8.1722939648034506e-20, 8.0517828972839211e-20, 7.9307883875099226e-20, | |
| /* [212] */ 7.8092777859524425e-20, 7.6872166028429042e-20, 7.5645683383965122e-20, 7.4412942930179128e-20, | |
| /* [216] */ 7.3173533545093332e-20, 7.1927017587631075e-20, 7.0672928197666785e-20, 6.9410766239500362e-20, | |
| /* [220] */ 6.8139996829256425e-20, 6.6860045374610234e-20, 6.5570293040210081e-20, 6.4270071533368528e-20, | |
| /* [224] */ 6.2958657080923559e-20, 6.1635263438143136e-20, 6.02990337321517e-20, 5.8949030892850181e-20, | |
| /* [228] */ 5.758422635988593e-20, 5.6203486669597397e-20, 5.4805557413499315e-20, 5.3389043909003295e-20, | |
| /* [232] */ 5.1952387717989917e-20, 5.0493837866338355e-20, 4.9011415222629489e-20, 4.7502867933366117e-20, | |
| /* [236] */ 4.5965615001265455e-20, 4.4396673897997565e-20, 4.2792566302148588e-20, 4.1149193273430015e-20, | |
| /* [240] */ 3.9461666762606287e-20, 3.7724077131401685e-20, 3.592916408620436e-20, 3.4067836691100565e-20, | |
| /* [244] */ 3.2128447641564046e-20, 3.0095646916399994e-20, 2.7948469455598328e-20, 2.5656913048718645e-20, | |
| /* [248] */ 2.3175209756803909e-20, 2.0426695228251291e-20, 1.7261770330213488e-20, 1.3281889259442579e-20, | |
| /* [252] */ 0, | |
| ] | |
| /// The precomputed ziggurat heights, denoted Y_i in the main text. | |
| private static let Y: [Double] = [ | |
| /* [ 0] */ 5.595205495112736e-23, 1.1802509982703313e-22, 1.8444423386735829e-22, 2.5439030466698309e-22, | |
| /* [ 4] */ 3.2737694311509334e-22, 4.0307732132706715e-22, 4.8125478319495115e-22, 5.6172914896583308e-22, | |
| /* [ 8] */ 6.4435820540443526e-22, 7.2902662343463681e-22, 8.1563888456321941e-22, 9.0411453683482223e-22, | |
| /* [ 12] */ 9.9438488486399206e-22, 1.0863906045969114e-21, 1.1800799775461269e-21, 1.2754075534831208e-21, | |
| /* [ 16] */ 1.372333117637729e-21, 1.4708208794375214e-21, 1.5708388257440445e-21, 1.6723581984374566e-21, | |
| /* [ 20] */ 1.7753530675030514e-21, 1.8797999785104595e-21, 1.9856776587832504e-21, 2.0929667704053244e-21, | |
| /* [ 24] */ 2.201649700995824e-21, 2.3117103852306179e-21, 2.4231341516125464e-21, 2.5359075901420891e-21, | |
| /* [ 28] */ 2.6500184374170538e-21, 2.7654554763660391e-21, 2.8822084483468604e-21, 3.0002679757547711e-21, | |
| /* [ 32] */ 3.1196254936130377e-21, 3.2402731888801749e-21, 3.3622039464187092e-21, 3.4854113007409036e-21, | |
| /* [ 36] */ 3.6098893927859475e-21, 3.7356329310971768e-21, 3.8626371568620053e-21, 3.9908978123552837e-21, | |
| /* [ 40] */ 4.1204111123918948e-21, 4.2511737184488913e-21, 4.3831827151633737e-21, 4.5164355889510656e-21, | |
| /* [ 44] */ 4.6509302085234806e-21, 4.7866648071096003e-21, 4.9236379662119969e-21, 5.0618486007478993e-21, | |
| /* [ 48] */ 5.2012959454434732e-21, 5.3419795423648946e-21, 5.4838992294830959e-21, 5.6270551301806347e-21, | |
| /* [ 52] */ 5.7714476436191935e-21, 5.9170774358950678e-21, 6.0639454319177027e-21, 6.2120528079531677e-21, | |
| /* [ 56] */ 6.3614009847804375e-21, 6.5119916214136427e-21, 6.6638266093481696e-21, 6.8169080672926277e-21, | |
| /* [ 60] */ 6.9712383363524377e-21, 7.1268199756340822e-21, 7.2836557582420336e-21, 7.4417486676430174e-21, | |
| /* [ 64] */ 7.6011018943746355e-21, 7.7617188330775411e-21, 7.9236030798322572e-21, 8.0867584297834842e-21, | |
| /* [ 68] */ 8.2511888750363333e-21, 8.4168986028103258e-21, 8.5838919938383098e-21, 8.7521736209986459e-21, | |
| /* [ 72] */ 8.9217482481700712e-21, 9.0926208292996504e-21, 9.2647965076751277e-21, 9.4382806153938292e-21, | |
| /* [ 76] */ 9.6130786730210328e-21, 9.7891963894314161e-21, 9.966639661827884e-21, 1.0145414575932636e-20, | |
| /* [ 80] */ 1.0325527406345955e-20, 1.0506984617068672e-20, 1.0689792862184811e-20, 1.0873958986701341e-20, | |
| /* [ 84] */ 1.10594900275424e-20, 1.1246393214695825e-20, 1.1434675972510121e-20, 1.1624345921140471e-20, | |
| /* [ 88] */ 1.1815410878142659e-20, 1.2007878860214202e-20, 1.2201758085082226e-20, 1.239705697353804e-20, | |
| /* [ 92] */ 1.2593784151618565e-20, 1.2791948452935152e-20, 1.29915589211506e-20, 1.3192624812605428e-20, | |
| /* [ 96] */ 1.3395155599094805e-20, 1.3599160970797774e-20, 1.3804650839360727e-20, 1.4011635341137284e-20, | |
| /* [100] */ 1.4220124840587164e-20, 1.4430129933836705e-20, 1.4641661452404201e-20, 1.485473046709328e-20, | |
| /* [104] */ 1.5069348292058084e-20, 1.5285526489044053e-20, 1.5503276871808626e-20, 1.5722611510726402e-20, | |
| /* [108] */ 1.5943542737583543e-20, 1.6166083150566702e-20, 1.6390245619451956e-20, 1.6616043290999594e-20, | |
| /* [112] */ 1.6843489594561079e-20, 1.7072598247904713e-20, 1.7303383263267072e-20, 1.7535858953637607e-20, | |
| /* [116] */ 1.7770039939284241e-20, 1.8005941154528286e-20, 1.8243577854777398e-20, 1.8482965623825808e-20, | |
| /* [120] */ 1.8724120381431627e-20, 1.8967058391181452e-20, 1.9211796268653192e-20, 1.9458350989888484e-20, | |
| /* [124] */ 1.9706739900186868e-20, 1.9956980723234356e-20, 2.0209091570579904e-20, 2.0463090951473895e-20, | |
| /* [128] */ 2.0718997783083593e-20, 2.097683140110135e-20, 2.123661157076213e-20, 2.1498358498287976e-20, | |
| /* [132] */ 2.1762092842777868e-20, 2.2027835728562592e-20, 2.2295608758045219e-20, 2.2565434025049041e-20, | |
| /* [136] */ 2.2837334128696004e-20, 2.311133218784001e-20, 2.3387451856080863e-20, 2.3665717337386111e-20, | |
| /* [140] */ 2.394615340234961e-20, 2.422878540511741e-20, 2.4513639301013211e-20, 2.4800741664897764e-20, | |
| /* [144] */ 2.5090119710298442e-20, 2.5381801309347597e-20, 2.56758150135705e-20, 2.5972190075566336e-20, | |
| /* [148] */ 2.6270956471628253e-20, 2.6572144925351523e-20, 2.6875786932281841e-20, 2.7181914785659148e-20, | |
| /* [152] */ 2.7490561603315974e-20, 2.7801761355793055e-20, 2.8115548895739172e-20, 2.8431959988666534e-20, | |
| /* [156] */ 2.8751031345137833e-20, 2.9072800654466307e-20, 2.9397306620015486e-20, 2.9724588996191657e-20, | |
| /* [160] */ 3.0054688627228112e-20, 3.0387647487867642e-20, 3.0723508726057078e-20, 3.1062316707775905e-20, | |
| /* [164] */ 3.1404117064129991e-20, 3.1748956740850969e-20, 3.2096884050352357e-20, 3.2447948726504914e-20, | |
| /* [168] */ 3.2802201982306013e-20, 3.3159696570631373e-20, 3.352048684827223e-20, 3.3884628843476888e-20, | |
| /* [172] */ 3.4252180327233346e-20, 3.4623200888548644e-20, 3.4997752014001677e-20, 3.537589717186906e-20, | |
| /* [176] */ 3.5757701901149035e-20, 3.6143233905835799e-20, 3.65325631548274e-20, 3.6925761987883572e-20, | |
| /* [180] */ 3.7322905228086981e-20, 3.7724070301302117e-20, 3.8129337363171041e-20, 3.8538789434235234e-20, | |
| /* [184] */ 3.8952512543827862e-20, 3.9370595883442399e-20, 3.9793131970351439e-20, 4.0220216822325769e-20, | |
| /* [188] */ 4.0651950144388133e-20, 4.1088435528630944e-20, 4.1529780668232712e-20, 4.1976097586926582e-20, | |
| /* [192] */ 4.2427502885307452e-20, 4.2884118005513604e-20, 4.3346069515987453e-20, 4.3813489418210257e-20, | |
| /* [196] */ 4.4286515477520838e-20, 4.4765291580372353e-20, 4.5249968120658306e-20, 4.5740702418054417e-20, | |
| /* [200] */ 4.6237659171683015e-20, 4.6741010952818368e-20, 4.7250938740823415e-20, 4.7767632507051219e-20, | |
| /* [204] */ 4.8291291852069895e-20, 4.8822126702292804e-20, 4.9360358072933852e-20, 4.9906218905182021e-20, | |
| /* [208] */ 5.0459954986625539e-20, 5.1021825965285324e-20, 5.1592106469178258e-20, 5.2171087345169234e-20, | |
| /* [212] */ 5.2759077033045284e-20, 5.3356403093325858e-20, 5.3963413910399511e-20, 5.4580480596259246e-20, | |
| /* [216] */ 5.5207999124535584e-20, 5.584639272987383e-20, 5.649611461419377e-20, 5.7157651009290713e-20, | |
| /* [220] */ 5.7831524654956632e-20, 5.8518298763794323e-20, 5.9218581558791713e-20, 5.99330314883387e-20, | |
| /* [224] */ 6.0662363246796887e-20, 6.1407354758435e-20, 6.2168855320499763e-20, 6.2947795150103727e-20, | |
| /* [228] */ 6.3745196643214394e-20, 6.4562187737537985e-20, 6.5400017881889097e-20, 6.6260077263309343e-20, | |
| /* [232] */ 6.714392014514662e-20, 6.8053293447301698e-20, 6.8990172088133e-20, 6.9956803158564498e-20, | |
| /* [236] */ 7.095576179487843e-20, 7.199002278894508e-20, 7.3063053739105458e-20, 7.4178938266266881e-20, | |
| /* [240] */ 7.5342542134173124e-20, 7.6559742171142969e-20, 7.783774986341285e-20, 7.9185582674029512e-20, | |
| /* [244] */ 8.06147755373533e-20, 8.2140502769818073e-20, 8.3783445978280519e-20, 8.5573129249678161e-20, | |
| /* [248] */ 8.75544596695901e-20, 8.9802388057706877e-20, 9.2462471421151086e-20, 9.5919641344951721e-20, | |
| /* [252] */ 1.0842021724855044e-19, | |
| ] | |
| /// Alias map for j in [0, 256). | |
| private static let MAP: [UInt8] = [ | |
| /* [ 0] */ 0, 0, 1, 235, 3, 4, 5, 0, | |
| /* [ 8] */ 0, 0, 0, 0, 0, 0, 0, 0, | |
| /* [ 16] */ 0, 0, 1, 1, 1, 1, 2, 2, | |
| /* [ 24] */ 252, 252, 252, 252, 252, 252, 252, 252, | |
| /* [ 32] */ 252, 252, 252, 252, 252, 252, 252, 252, | |
| /* [ 40] */ 252, 252, 252, 252, 252, 252, 252, 252, | |
| /* [ 48] */ 252, 252, 252, 252, 252, 252, 252, 252, | |
| /* [ 56] */ 252, 252, 252, 252, 252, 252, 252, 252, | |
| /* [ 64] */ 252, 252, 252, 252, 252, 252, 252, 252, | |
| /* [ 72] */ 252, 252, 252, 252, 252, 252, 252, 252, | |
| /* [ 80] */ 252, 252, 252, 252, 252, 252, 252, 252, | |
| /* [ 88] */ 252, 252, 252, 252, 252, 252, 252, 252, | |
| /* [ 96] */ 252, 252, 252, 252, 252, 252, 252, 252, | |
| /* [104] */ 252, 252, 252, 252, 252, 252, 252, 252, | |
| /* [112] */ 252, 252, 252, 252, 252, 252, 252, 252, | |
| /* [120] */ 252, 252, 252, 252, 252, 252, 252, 252, | |
| /* [128] */ 252, 252, 252, 252, 252, 252, 252, 252, | |
| /* [136] */ 252, 252, 252, 252, 252, 252, 252, 252, | |
| /* [144] */ 252, 252, 252, 252, 252, 252, 252, 252, | |
| /* [152] */ 252, 252, 252, 252, 252, 252, 252, 252, | |
| /* [160] */ 252, 252, 252, 252, 252, 252, 252, 252, | |
| /* [168] */ 252, 252, 252, 252, 252, 252, 252, 252, | |
| /* [176] */ 252, 251, 251, 251, 251, 251, 251, 251, | |
| /* [184] */ 251, 251, 251, 251, 251, 251, 250, 250, | |
| /* [192] */ 250, 250, 250, 250, 250, 249, 249, 249, | |
| /* [200] */ 249, 249, 249, 248, 248, 248, 248, 247, | |
| /* [208] */ 247, 247, 247, 246, 246, 246, 245, 245, | |
| /* [216] */ 244, 244, 243, 243, 242, 241, 241, 240, | |
| /* [224] */ 239, 237, 3, 3, 4, 4, 6, 0, | |
| /* [232] */ 0, 0, 0, 236, 237, 238, 239, 240, | |
| /* [240] */ 241, 242, 243, 244, 245, 246, 247, 248, | |
| /* [248] */ 249, 250, 251, 252, 2, 0, 0, 0, | |
| ] | |
| /// The alias inverse PMF for j in-place of j. | |
| private static let IPMF: [Int64] = [ | |
| /* [ 0] */ 9223372036854774016, 1623796909450834944, 2664290944894291200, 7387971354164060928, | |
| /* [ 4] */ 6515064486552723200, 8840508362680718848, 6099647593382936320, 7673130333659513856, | |
| /* [ 8] */ 6220332867583438080, 5045979640552813824, 4075305837223955456, 3258413672162525440, | |
| /* [ 12] */ 2560664887087762432, 1957224924672899584, 1429800935350577408, 964606309710808320, | |
| /* [ 16] */ 551043923599587072, 180827629096890368, -152619738120023552, -454588624410291456, | |
| /* [ 20] */ -729385126147774976, -980551509819447040, -1211029700667463936, -1423284293868548352, | |
| /* [ 24] */ -1619396356369050368, -1801135830956211712, -1970018048575618048, -2127348289059705344, | |
| /* [ 28] */ -2274257249303686400, -2411729520096655360, -2540626634159181056, -2661705860113406464, | |
| /* [ 32] */ -2775635634532450560, -2883008316030465280, -2984350790383654912, -3080133339198116352, | |
| /* [ 36] */ -3170777096303091200, -3256660348483819008, -3338123885075136256, -3415475560473299200, | |
| /* [ 40] */ -3488994201966428160, -3558932970354473216, -3625522261068041216, -3688972217741989376, | |
| /* [ 44] */ -3749474917563782656, -3807206277531056128, -3862327722496843520, -3914987649156779776, | |
| /* [ 48] */ -3965322714631865344, -4013458973776895488, -4059512885612783360, -4103592206186241024, | |
| /* [ 52] */ -4145796782586128128, -4186219260694347008, -4224945717447275264, -4262056226866285568, | |
| /* [ 56] */ -4297625367836519680, -4331722680528537344, -4364413077437472512, -4395757214229401600, | |
| /* [ 60] */ -4425811824915135744, -4454630025296932608, -4482261588141290496, -4508753193105288192, | |
| /* [ 64] */ -4534148654077808896, -4558489126279958272, -4581813295192216576, -4604157549138257664, | |
| /* [ 68] */ -4625556137145255168, -4646041313519104512, -4665643470413305856, -4684391259530326528, | |
| /* [ 72] */ -4702311703971761664, -4719430301145103360, -4735771117539946240, -4751356876102087168, | |
| /* [ 76] */ -4766209036859133952, -4780347871386013440, -4793792531638892032, -4806561113635132672, | |
| /* [ 80] */ -4818670716409306624, -4830137496634465536, -4840976719260837888, -4851202804490348800, | |
| /* [ 84] */ -4860829371376460032, -4869869278311657472, -4878334660640771072, -4886236965617427200, | |
| /* [ 88] */ -4893586984900802560, -4900394884772702720, -4906670234238885376, -4912422031164496896, | |
| /* [ 92] */ -4917658726580119808, -4922388247283532288, -4926618016851066624, -4930354975163335168, | |
| /* [ 96] */ -4933605596540651264, -4936375906575303936, -4938671497741366016, -4940497543854575616, | |
| /* [100] */ -4941858813449629440, -4942759682136114944, -4943204143989086720, -4943195822025528064, | |
| /* [104] */ -4942737977813206528, -4941833520255033344, -4940485013586738944, -4938694684624359424, | |
| /* [108] */ -4936464429291795968, -4933795818458825728, -4930690103114057984, -4927148218896864000, | |
| /* [112] */ -4923170790008275968, -4918758132519213568, -4913910257091645696, -4908626871126539264, | |
| /* [116] */ -4902907380349533952, -4896750889844272896, -4890156204540531200, -4883121829162554368, | |
| /* [120] */ -4875645967641781248, -4867726521994927104, -4859361090668103424, -4850546966345113600, | |
| /* [124] */ -4841281133215539200, -4831560263698491904, -4821380714613447424, -4810738522790066176, | |
| /* [128] */ -4799629400105481984, -4788048727936307200, -4775991551010514944, -4763452570642114304, | |
| /* [132] */ -4750426137329494528, -4736906242696389120, -4722886510751377664, -4708360188440089088, | |
| /* [136] */ -4693320135461421056, -4677758813316108032, -4661668273553489152, -4645040145179241472, | |
| /* [140] */ -4627865621182772224, -4610135444140930048, -4591839890849345536, -4572968755929961472, | |
| /* [144] */ -4553511334358205696, -4533456402849101568, -4512792200036279040, -4491506405372580864, | |
| /* [148] */ -4469586116675402496, -4447017826233107968, -4423787395382284800, -4399880027458416384, | |
| /* [152] */ -4375280239014115072, -4349971829190472192, -4323937847117721856, -4297160557210933504, | |
| /* [156] */ -4269621402214949888, -4241300963840749312, -4212178920821861632, -4182234004204451584, | |
| /* [160] */ -4151443949668877312, -4119785446662287616, -4087234084103201536, -4053764292396156928, | |
| /* [164] */ -4019349281473081856, -3983960974549692672, -3947569937258423296, -3910145301787345664, | |
| /* [168] */ -3871654685619032064, -3832064104425388800, -3791337878631544832, -3749438533114327552, | |
| /* [172] */ -3706326689447984384, -3661960950051848192, -3616297773528534784, -3569291340409189376, | |
| /* [176] */ -3520893408440946176, -3471053156460654336, -3419717015797782528, -3366828488034805504, | |
| /* [180] */ -3312327947826460416, -3256152429334010368, -3198235394669719040, -3138506482563172864, | |
| /* [184] */ -3076891235255162880, -3013310801389730816, -2947681612411374848, -2879915029671670784, | |
| /* [188] */ -2809916959107513856, -2737587429961866240, -2662820133571325696, -2585501917733380096, | |
| /* [192] */ -2505512231579385344, -2422722515205211648, -2336995527534088448, -2248184604988727552, | |
| /* [196] */ -2156132842510765056, -2060672187261025536, -1961622433929371904, -1858790108950105600, | |
| /* [200] */ -1751967229002895616, -1640929916937142784, -1525436855617582592, -1405227557075253248, | |
| /* [204] */ -1280020420662650112, -1149510549536596224, -1013367289578704896, -871231448632104192, | |
| /* [208] */ -722712146453667840, -567383236774436096, -404779231966938368, -234390647591545856, | |
| /* [212] */ -55658667960119296, 132030985907841280, 329355128892811776, 537061298001085184, | |
| /* [216] */ 755977262693564160, 987022116608033280, 1231219266829431296, 1489711711346518528, | |
| /* [220] */ 1763780090187553792, 2054864117341795072, 2364588157623768832, 2694791916990503168, | |
| /* [224] */ 3047567482883476224, 3425304305830816256, 3830744187097297920, 4267048975685830400, | |
| /* [228] */ 4737884547990017280, 5247525842198998272, 5800989391535355392, 6404202162993295360, | |
| /* [232] */ 7064218894258540544, 7789505049452331520, 8590309807749444864, 7643763810684489984, | |
| /* [236] */ 8891950541491446016, 5457384281016206080, 9083704440929284096, 7976211653914433280, | |
| /* [240] */ 8178631350487117568, 2821287825726744832, 6322989683301709568, 4309503753387611392, | |
| /* [244] */ 4685170734960170496, 8404845967535199744, 7330522972447554048, 1960945799076992000, | |
| /* [248] */ 4742910674644899072, -751799822533509888, 7023456603741959936, 3843116882594676224, | |
| /* [252] */ 3927231442413903104, -9223372036854775808, -9223372036854775808, -9223372036854775808, | |
| ] | |
| var rng: UniformRandomProvider | |
| init(rng: UniformRandomProvider) { | |
| self.rng = rng | |
| } | |
| mutating func sample() -> Double { | |
| let x = rng.nextUniform() | |
| let i = Int(truncatingIfNeeded: x) & Self.maskInt8 | |
| if i < Self.I_MAX { | |
| let magnitude = Double(x >> 1) | |
| return Self.X[i] * magnitude | |
| } | |
| return edgeSample(initial: x) | |
| } | |
| private mutating func edgeSample(initial xx: UInt64) -> Double { | |
| var j = selectRegion() | |
| if j != 0 { | |
| return sampleOverhang(j: j, initial: xx) | |
| } | |
| var x0 = Self.X_0 | |
| while true { | |
| let x = rng.nextUniform() | |
| let i = Int(truncatingIfNeeded: x) & Self.maskInt8 | |
| if i < Self.I_MAX { | |
| let magnitude = Double(x >> 1) | |
| return x0 + Self.X[i] * magnitude | |
| } | |
| j = selectRegion() | |
| if j != 0 { | |
| return x0 + sampleOverhang(j: j, initial: x) | |
| } | |
| x0 += Self.X_0 | |
| } | |
| } | |
| private mutating func selectRegion() -> Int { | |
| let x = rng.nextUniform() | |
| let j = Int(truncatingIfNeeded: x) & Self.maskInt8 | |
| let alias = Int(Self.MAP[j]) & Self.maskInt8 | |
| let signedX = Int64(bitPattern: x) | |
| return signedX >= Self.IPMF[j] ? alias : j | |
| } | |
| private mutating func sampleOverhang(j: Int, initial xx: UInt64) -> Double { | |
| var u1 = xx >> 1 | |
| while true { | |
| let u = randomInt63() | |
| let u2: UInt64 | |
| if u1 < u { | |
| u2 = u | |
| } else { | |
| u2 = u1 | |
| u1 = u | |
| } | |
| let x = Self.interpolate(Self.X, j: j, u: u1) | |
| if u2 &- u1 >= Self.E_MAX { | |
| return x | |
| } | |
| if Self.interpolate(Self.Y, j: j, u: u2) <= exp(-x) { | |
| return x | |
| } | |
| u1 = randomInt63() | |
| } | |
| } | |
| } | |
| } | |
| extension ZigguratSamplers { | |
| /// Modified ziggurat method for sampling from a Gaussian distribution with mean 0 and standard deviation 1. | |
| struct NormalizedGaussian: ZigguratSampler { | |
| private static let I_MAX = 253 | |
| private static let J_INFLECTION = 204 | |
| private static let CONVEX_E_MAX: Int64 = -2269182951627976004 | |
| private static let CONCAVE_E_MAX: UInt64 = 760463704284035184 | |
| private static let X_0 = 3.6360066255009455861 | |
| private static let ONE_OVER_X_0 = 1.0 / X_0 | |
| private static let X: [Double] = [ | |
| /* [ 0] */ 3.9421662825398133e-19, 3.7204945004119012e-19, 3.5827024480628678e-19, 3.4807476236540249e-19, | |
| /* [ 4] */ 3.3990177171882136e-19, 3.3303778360340139e-19, 3.270943881761755e-19, 3.21835771324951e-19, | |
| /* [ 8] */ 3.1710758541840432e-19, 3.1280307407034065e-19, 3.0884520655804019e-19, 3.0517650624107352e-19, | |
| /* [ 12] */ 3.01752902925846e-19, 2.985398344070532e-19, 2.9550967462801797e-19, 2.9263997988491663e-19, | |
| /* [ 16] */ 2.8991225869977476e-19, 2.8731108780226291e-19, 2.8482346327101335e-19, 2.8243831535194389e-19, | |
| /* [ 20] */ 2.8014613964727031e-19, 2.7793871261807797e-19, 2.7580886921411212e-19, 2.7375032698308758e-19, | |
| /* [ 24] */ 2.7175754543391047e-19, 2.6982561247538484e-19, 2.6795015188771505e-19, 2.6612724730440033e-19, | |
| /* [ 28] */ 2.6435337927976633e-19, 2.6262537282028438e-19, 2.6094035335224142e-19, 2.5929570954331002e-19, | |
| /* [ 32] */ 2.5768906173214726e-19, 2.5611823497719608e-19, 2.5458123593393361e-19, 2.5307623292372459e-19, | |
| /* [ 36] */ 2.51601538677984e-19, 2.5015559533646191e-19, 2.4873696135403158e-19, 2.4734430003079206e-19, | |
| /* [ 40] */ 2.4597636942892726e-19, 2.446320134791245e-19, 2.4331015411139206e-19, 2.4200978427132955e-19, | |
| /* [ 44] */ 2.4072996170445879e-19, 2.3946980340903347e-19, 2.3822848067252674e-19, 2.3700521461931801e-19, | |
| /* [ 48] */ 2.357992722074133e-19, 2.3460996262069972e-19, 2.3343663401054455e-19, 2.322786705467384e-19, | |
| /* [ 52] */ 2.3113548974303765e-19, 2.3000654002704238e-19, 2.2889129852797606e-19, 2.2778926905921897e-19, | |
| /* [ 56] */ 2.2669998027527321e-19, 2.2562298398527416e-19, 2.245578536072726e-19, 2.2350418274933911e-19, | |
| /* [ 60] */ 2.2246158390513294e-19, 2.2142968725296249e-19, 2.2040813954857555e-19, 2.1939660310297601e-19, | |
| /* [ 64] */ 2.1839475483749618e-19, 2.1740228540916853e-19, 2.1641889840016519e-19, 2.1544430956570613e-19, | |
| /* [ 68] */ 2.1447824613540345e-19, 2.1352044616350571e-19, 2.1257065792395107e-19, 2.1162863934653125e-19, | |
| /* [ 72] */ 2.1069415749082026e-19, 2.0976698805483467e-19, 2.0884691491567363e-19, 2.0793372969963634e-19, | |
| /* [ 76] */ 2.0702723137954107e-19, 2.0612722589717129e-19, 2.0523352580895635e-19, 2.0434594995315797e-19, | |
| /* [ 80] */ 2.0346432313698148e-19, 2.0258847584216418e-19, 2.0171824394771313e-19, 2.0085346846857531e-19, | |
| /* [ 84] */ 1.9999399530912015e-19, 1.9913967503040585e-19, 1.9829036263028144e-19, 1.9744591733545175e-19, | |
| /* [ 88] */ 1.9660620240469857e-19, 1.9577108494251485e-19, 1.9494043572246307e-19, 1.9411412901962161e-19, | |
| /* [ 92] */ 1.9329204245152935e-19, 1.9247405682708168e-19, 1.9166005600287074e-19, 1.9084992674649826e-19, | |
| /* [ 96] */ 1.900435586064234e-19, 1.8924084378793725e-19, 1.8844167703488436e-19, 1.8764595551677749e-19, | |
| /* [100] */ 1.868535787209745e-19, 1.8606444834960934e-19, 1.8527846822098793e-19, 1.8449554417517928e-19, | |
| /* [104] */ 1.8371558398354868e-19, 1.8293849726199566e-19, 1.8216419538767393e-19, 1.8139259141898448e-19, | |
| /* [108] */ 1.8062360001864453e-19, 1.7985713737964743e-19, 1.7909312115393845e-19, 1.78331470383642e-19, | |
| /* [112] */ 1.7757210543468428e-19, 1.7681494793266395e-19, 1.760599207008314e-19, 1.7530694770004409e-19, | |
| /* [116] */ 1.7455595397057217e-19, 1.7380686557563475e-19, 1.7305960954655264e-19, 1.7231411382940904e-19, | |
| /* [120] */ 1.7157030723311378e-19, 1.7082811937877138e-19, 1.7008748065025788e-19, 1.6934832214591352e-19, | |
| /* [124] */ 1.6861057563126349e-19, 1.6787417349268046e-19, 1.6713904869190636e-19, 1.6640513472135291e-19, | |
| /* [128] */ 1.6567236556010242e-19, 1.6494067563053266e-19, 1.6420999975549115e-19, 1.6348027311594532e-19, | |
| /* [132] */ 1.6275143120903661e-19, 1.6202340980646725e-19, 1.6129614491314931e-19, 1.6056957272604589e-19, | |
| /* [136] */ 1.5984362959313479e-19, 1.5911825197242491e-19, 1.5839337639095554e-19, 1.57668939403708e-19, | |
| /* [140] */ 1.5694487755235889e-19, 1.5622112732380261e-19, 1.554976251083707e-19, 1.5477430715767271e-19, | |
| /* [144] */ 1.540511095419833e-19, 1.5332796810709688e-19, 1.5260481843056974e-19, 1.5188159577726683e-19, | |
| /* [148] */ 1.5115823505412761e-19, 1.5043467076406199e-19, 1.4971083695888395e-19, 1.4898666719118714e-19, | |
| /* [152] */ 1.4826209446506113e-19, 1.4753705118554365e-19, 1.468114691066983e-19, 1.4608527927820112e-19, | |
| /* [156] */ 1.4535841199031451e-19, 1.4463079671711862e-19, 1.4390236205786415e-19, 1.4317303567630177e-19, | |
| /* [160] */ 1.4244274423783481e-19, 1.4171141334433217e-19, 1.4097896746642792e-19, 1.4024532987312287e-19, | |
| /* [164] */ 1.3951042255849034e-19, 1.3877416616527576e-19, 1.3803647990516385e-19, 1.3729728147547174e-19, | |
| /* [168] */ 1.3655648697200824e-19, 1.3581401079782068e-19, 1.3506976556752901e-19, 1.3432366200692418e-19, | |
| /* [172] */ 1.3357560884748263e-19, 1.3282551271542047e-19, 1.3207327801488087e-19, 1.3131880680481524e-19, | |
| /* [176] */ 1.3056199866908076e-19, 1.2980275057923788e-19, 1.2904095674948608e-19, 1.2827650848312727e-19, | |
| /* [180] */ 1.2750929400989213e-19, 1.2673919831340482e-19, 1.2596610294799512e-19, 1.2518988584399374e-19, | |
| /* [184] */ 1.2441042110056523e-19, 1.2362757876504165e-19, 1.2284122459762072e-19, 1.2205121982017852e-19, | |
| /* [188] */ 1.2125742084782245e-19, 1.2045967900166973e-19, 1.196578402011802e-19, 1.1885174463419555e-19, | |
| /* [192] */ 1.1804122640264091e-19, 1.1722611314162064e-19, 1.1640622560939109e-19, 1.1558137724540874e-19, | |
| /* [196] */ 1.1475137369333185e-19, 1.1391601228549047e-19, 1.1307508148492592e-19, 1.1222836028063025e-19, | |
| /* [200] */ 1.1137561753107903e-19, 1.1051661125053526e-19, 1.0965108783189755e-19, 1.0877878119905372e-19, | |
| /* [204] */ 1.0789941188076655e-19, 1.070126859970364e-19, 1.0611829414763286e-19, 1.0521591019102928e-19, | |
| /* [208] */ 1.0430518990027552e-19, 1.0338576948035472e-19, 1.0245726392923699e-19, 1.015192652220931e-19, | |
| /* [212] */ 1.0057134029488235e-19, 9.9613028799672809e-20, 9.8643840599459914e-20, 9.7663252964755816e-20, | |
| /* [216] */ 9.6670707427623454e-20, 9.566560624086667e-20, 9.4647308380433213e-20, 9.3615125017323508e-20, | |
| /* [220] */ 9.2568314370887282e-20, 9.1506075837638774e-20, 9.0427543267725716e-20, 8.933177723376368e-20, | |
| /* [224] */ 8.8217756102327883e-20, 8.7084365674892319e-20, 8.5930387109612162e-20, 8.4754482764244349e-20, | |
| /* [228] */ 8.3555179508462343e-20, 8.2330848933585364e-20, 8.1079683729129853e-20, 7.9799669284133864e-20, | |
| /* [232] */ 7.8488549286072745e-20, 7.7143783700934692e-20, 7.5762496979467566e-20, 7.4341413578485329e-20, | |
| /* [236] */ 7.2876776807378431e-20, 7.1364245443525374e-20, 6.9798760240761066e-20, 6.8174368944799054e-20, | |
| /* [240] */ 6.6483992986198539e-20, 6.4719110345162767e-20, 6.2869314813103699e-20, 6.0921687548281263e-20, | |
| /* [244] */ 5.8859873575576818e-20, 5.6662675116090981e-20, 5.4301813630894571e-20, 5.173817174449422e-20, | |
| /* [248] */ 4.8915031722398545e-20, 4.5744741890755301e-20, 4.2078802568583416e-20, 3.7625986722404761e-20, | |
| /* [252] */ 3.1628589805881879e-20, 0, | |
| ] | |
| private static let Y: [Double] = [ | |
| /* [ 0] */ 1.4598410796619063e-22, 3.0066613427942797e-22, 4.6129728815103466e-22, 6.2663350049234362e-22, | |
| /* [ 4] */ 7.9594524761881544e-22, 9.6874655021705039e-22, 1.1446877002379439e-21, 1.3235036304379167e-21, | |
| /* [ 8] */ 1.5049857692053131e-21, 1.6889653000719298e-21, 1.8753025382711626e-21, 2.0638798423695191e-21, | |
| /* [ 12] */ 2.2545966913644708e-21, 2.4473661518801799e-21, 2.6421122727763533e-21, 2.8387681187879908e-21, | |
| /* [ 16] */ 3.0372742567457284e-21, 3.2375775699986589e-21, 3.439630315794878e-21, 3.6433893657997798e-21, | |
| /* [ 20] */ 3.8488155868912312e-21, 4.0558733309492775e-21, 4.264530010428359e-21, 4.4747557422305067e-21, | |
| /* [ 24] */ 4.6865230465355582e-21, 4.8998065902775257e-21, 5.1145829672105489e-21, 5.3308305082046173e-21, | |
| /* [ 28] */ 5.5485291167031758e-21, 5.7676601252690476e-21, 5.9882061699178461e-21, 6.2101510795442221e-21, | |
| /* [ 32] */ 6.4334797782257209e-21, 6.6581781985713897e-21, 6.8842332045893181e-21, 7.1116325227957095e-21, | |
| /* [ 36] */ 7.3403646804903092e-21, 7.5704189502886418e-21, 7.8017853001379744e-21, 8.0344543481570017e-21, | |
| /* [ 40] */ 8.2684173217333118e-21, 8.5036660203915022e-21, 8.7401927820109521e-21, 8.9779904520281901e-21, | |
| /* [ 44] */ 9.2170523553061439e-21, 9.457372270392882e-21, 9.698944405926943e-21, 9.9417633789758424e-21, | |
| /* [ 48] */ 1.0185824195119818e-20, 1.043112223011477e-20, 1.0677653212987396e-20, 1.0925413210432004e-20, | |
| /* [ 52] */ 1.1174398612392891e-20, 1.1424606118728715e-20, 1.1676032726866302e-20, 1.1928675720361027e-20, | |
| /* [ 56] */ 1.2182532658289373e-20, 1.2437601365406785e-20, 1.2693879923010674e-20, 1.2951366660454145e-20, | |
| /* [ 60] */ 1.3210060147261461e-20, 1.3469959185800733e-20, 1.3731062804473644e-20, 1.3993370251385596e-20, | |
| /* [ 64] */ 1.4256880988463136e-20, 1.4521594685988369e-20, 1.4787511217522902e-20, 1.505463065519617e-20, | |
| /* [ 68] */ 1.5322953265335218e-20, 1.5592479504415048e-20, 1.5863210015310328e-20, 1.6135145623830982e-20, | |
| /* [ 72] */ 1.6408287335525592e-20, 1.6682636332737932e-20, 1.6958193971903124e-20, 1.7234961781071113e-20, | |
| /* [ 76] */ 1.7512941457646084e-20, 1.7792134866331487e-20, 1.807254403727107e-20, 1.8354171164377277e-20, | |
| /* [ 80] */ 1.8637018603838945e-20, 1.8921088872801004e-20, 1.9206384648209468e-20, 1.9492908765815636e-20, | |
| /* [ 84] */ 1.9780664219333857e-20, 2.0069654159747839e-20, 2.0359881894760859e-20, 2.0651350888385696e-20, | |
| /* [ 88] */ 2.0944064760670539e-20, 2.1238027287557466e-20, 2.1533242400870487e-20, 2.1829714188430474e-20, | |
| /* [ 92] */ 2.2127446894294597e-20, 2.242644491911827e-20, 2.2726712820637798e-20, 2.3028255314272276e-20, | |
| /* [ 96] */ 2.3331077273843558e-20, 2.3635183732413286e-20, 2.3940579883236352e-20, 2.4247271080830277e-20, | |
| /* [100] */ 2.455526284216033e-20, 2.4864560847940368e-20, 2.5175170944049622e-20, 2.5487099143065929e-20, | |
| /* [104] */ 2.5800351625915997e-20, 2.6114934743643687e-20, 2.6430855019297323e-20, 2.6748119149937411e-20, | |
| /* [108] */ 2.7066734008766247e-20, 2.7386706647381193e-20, 2.7708044298153558e-20, 2.8030754376735269e-20, | |
| /* [112] */ 2.8354844484695747e-20, 2.8680322412291631e-20, 2.9007196141372126e-20, 2.9335473848423219e-20, | |
| /* [116] */ 2.9665163907753988e-20, 2.9996274894828624e-20, 3.0328815589748056e-20, 3.0662794980885287e-20, | |
| /* [120] */ 3.099822226867876e-20, 3.1335106869588609e-20, 3.1673458420220558e-20, 3.2013286781622988e-20, | |
| /* [124] */ 3.2354602043762612e-20, 3.2697414530184806e-20, 3.304173480286495e-20, 3.3387573667257349e-20, | |
| /* [128] */ 3.3734942177548938e-20, 3.4083851642125208e-20, 3.4434313629256243e-20, 3.4786339973011376e-20, | |
| /* [132] */ 3.5139942779411164e-20, 3.5495134432826171e-20, 3.585192760263246e-20, 3.6210335250134172e-20, | |
| /* [136] */ 3.6570370635764384e-20, 3.6932047326575882e-20, 3.7295379204034252e-20, 3.7660380472126401e-20, | |
| /* [140] */ 3.8027065665798284e-20, 3.8395449659736649e-20, 3.8765547677510167e-20, 3.9137375301086406e-20, | |
| /* [144] */ 3.9510948480742172e-20, 3.988628354538543e-20, 4.0263397213308566e-20, 4.0642306603393541e-20, | |
| /* [148] */ 4.1023029246790967e-20, 4.1405583099096438e-20, 4.1789986553048817e-20, 4.2176258451776819e-20, | |
| /* [152] */ 4.2564418102621759e-20, 4.2954485291566197e-20, 4.3346480298300118e-20, 4.3740423911958146e-20, | |
| /* [156] */ 4.4136337447563716e-20, 4.4534242763218286e-20, 4.4934162278076256e-20, 4.5336118991149025e-20, | |
| /* [160] */ 4.5740136500984466e-20, 4.6146239026271279e-20, 4.6554451427421133e-20, 4.6964799229185088e-20, | |
| /* [164] */ 4.7377308644364938e-20, 4.7792006598684169e-20, 4.8208920756888113e-20, 4.8628079550147814e-20, | |
| /* [168] */ 4.9049512204847653e-20, 4.9473248772842596e-20, 4.9899320163277674e-20, 5.0327758176068971e-20, | |
| /* [172] */ 5.0758595537153414e-20, 5.1191865935622696e-20, 5.1627604062866059e-20, 5.2065845653856416e-20, | |
| /* [176] */ 5.2506627530725194e-20, 5.2949987648783448e-20, 5.3395965145159426e-20, 5.3844600390237576e-20, | |
| /* [180] */ 5.4295935042099358e-20, 5.4750012104183868e-20, 5.5206875986405073e-20, 5.5666572569983821e-20, | |
| /* [184] */ 5.6129149276275792e-20, 5.6594655139902476e-20, 5.7063140886520563e-20, 5.7534659015596918e-20, | |
| /* [188] */ 5.8009263888591218e-20, 5.8487011822987583e-20, 5.8967961192659803e-20, 5.9452172535103471e-20, | |
| /* [192] */ 5.9939708666122605e-20, 6.0430634802618929e-20, 6.0925018694200531e-20, 6.142293076440286e-20, | |
| /* [196] */ 6.1924444262401531e-20, 6.2429635426193939e-20, 6.2938583658336214e-20, 6.3451371715447563e-20, | |
| /* [200] */ 6.3968085912834963e-20, 6.4488816345752736e-20, 6.5013657128995346e-20, 6.5542706656731714e-20, | |
| /* [204] */ 6.6076067884730717e-20, 6.6613848637404196e-20, 6.715616194241298e-20, 6.770312639595058e-20, | |
| /* [208] */ 6.8254866562246408e-20, 6.8811513411327825e-20, 6.9373204799659681e-20, 6.9940085998959109e-20, | |
| /* [212] */ 7.0512310279279503e-20, 7.1090039553397167e-20, 7.1673445090644796e-20, 7.2262708309655784e-20, | |
| /* [216] */ 7.2858021661057338e-20, 7.34595896130358e-20, 7.4067629754967553e-20, 7.4682374037052817e-20, | |
| /* [220] */ 7.5304070167226666e-20, 7.5932983190698547e-20, 7.6569397282483754e-20, 7.7213617789487678e-20, | |
| /* [224] */ 7.7865973566417016e-20, 7.8526819659456755e-20, 7.919654040385056e-20, 7.9875553017037968e-20, | |
| /* [228] */ 8.056431178890163e-20, 8.1263312996426176e-20, 8.1973100703706304e-20, 8.2694273652634034e-20, | |
| /* [232] */ 8.3427493508836792e-20, 8.4173494807453416e-20, 8.4933097052832066e-20, 8.5707219578230905e-20, | |
| /* [236] */ 8.6496899985930695e-20, 8.7303317295655327e-20, 8.8127821378859504e-20, 8.8971970928196666e-20, | |
| /* [240] */ 8.9837583239314064e-20, 9.0726800697869543e-20, 9.1642181484063544e-20, 9.2586826406702765e-20, | |
| /* [244] */ 9.3564561480278864e-20, 9.4580210012636175e-20, 9.5640015550850358e-20, 9.675233477050313e-20, | |
| /* [248] */ 9.7928851697808831e-20, 9.9186905857531331e-20, 1.0055456271343397e-19, 1.0208407377305566e-19, | |
| /* [252] */ 1.0390360993240711e-19, 1.0842021724855044e-19, | |
| ] | |
| private static let MAP: [UInt8] = [ | |
| /* [ 0] */ 0, 0, 239, 2, 0, 0, 0, 0, | |
| /* [ 8] */ 0, 0, 0, 0, 1, 1, 1, 253, | |
| /* [ 16] */ 253, 253, 253, 253, 253, 253, 253, 253, | |
| /* [ 24] */ 253, 253, 253, 253, 253, 253, 253, 253, | |
| /* [ 32] */ 253, 253, 253, 253, 253, 253, 253, 253, | |
| /* [ 40] */ 253, 253, 253, 253, 253, 253, 253, 253, | |
| /* [ 48] */ 253, 253, 253, 253, 253, 253, 253, 253, | |
| /* [ 56] */ 253, 253, 253, 253, 253, 253, 253, 253, | |
| /* [ 64] */ 253, 253, 253, 253, 253, 253, 253, 253, | |
| /* [ 72] */ 253, 253, 253, 253, 253, 253, 253, 253, | |
| /* [ 80] */ 253, 253, 253, 253, 253, 253, 253, 253, | |
| /* [ 88] */ 253, 253, 253, 253, 253, 253, 253, 253, | |
| /* [ 96] */ 253, 253, 253, 253, 253, 253, 253, 253, | |
| /* [104] */ 253, 253, 253, 253, 253, 253, 253, 253, | |
| /* [112] */ 253, 253, 253, 253, 253, 253, 253, 253, | |
| /* [120] */ 253, 253, 253, 253, 253, 253, 253, 253, | |
| /* [128] */ 253, 253, 253, 253, 253, 253, 253, 253, | |
| /* [136] */ 253, 253, 253, 253, 253, 253, 253, 253, | |
| /* [144] */ 253, 253, 253, 253, 253, 253, 253, 253, | |
| /* [152] */ 253, 253, 253, 253, 253, 253, 253, 253, | |
| /* [160] */ 253, 253, 253, 253, 253, 253, 253, 253, | |
| /* [168] */ 253, 253, 253, 253, 253, 253, 253, 253, | |
| /* [176] */ 253, 253, 253, 253, 253, 253, 253, 253, | |
| /* [184] */ 253, 253, 253, 253, 253, 253, 253, 253, | |
| /* [192] */ 253, 253, 252, 252, 252, 252, 252, 252, | |
| /* [200] */ 252, 252, 252, 252, 251, 251, 251, 251, | |
| /* [208] */ 251, 251, 251, 250, 250, 250, 250, 250, | |
| /* [216] */ 249, 249, 249, 248, 248, 248, 247, 247, | |
| /* [224] */ 247, 246, 246, 245, 244, 244, 243, 242, | |
| /* [232] */ 240, 2, 2, 3, 3, 0, 0, 240, | |
| /* [240] */ 241, 242, 243, 244, 245, 246, 247, 248, | |
| /* [248] */ 249, 250, 251, 252, 253, 1, 0, 0, | |
| ] | |
| private static let IPMF: [Int64] = [ | |
| /* [ 0] */ 9223372036854775296, 1100243796534090752, 7866600928998383104, 6788754710675124736, | |
| /* [ 4] */ 9022865200181688320, 6522434035205502208, 4723064097360024576, 3360495653216416000, | |
| /* [ 8] */ 2289663232373870848, 1423968905551920384, 708364817827798016, 106102487305601280, | |
| /* [ 12] */ -408333464665794560, -853239722779025152, -1242095211825521408, -1585059631105762048, | |
| /* [ 16] */ -1889943050287169024, -2162852901990669824, -2408637386594511104, -2631196530262954496, | |
| /* [ 20] */ -2833704942520925696, -3018774289025787392, -3188573753472222208, -3344920681707410944, | |
| /* [ 24] */ -3489349705062150656, -3623166100042179584, -3747487436868335360, -3863276422712173824, | |
| /* [ 28] */ -3971367044063130880, -4072485557029824000, -4167267476830916608, -4256271432240159744, | |
| /* [ 32] */ -4339990541927306752, -4418861817133802240, -4493273980372377088, -4563574004462246656, | |
| /* [ 36] */ -4630072609770453760, -4693048910430964992, -4752754358862894848, -4809416110052769536, | |
| /* [ 40] */ -4863239903586985984, -4914412541515875840, -4963104028439161088, -5009469424769119232, | |
| /* [ 44] */ -5053650458856559360, -5095776932695077632, -5135967952544929024, -5174333008451230720, | |
| /* [ 48] */ -5210972924952654336, -5245980700100460288, -5279442247516297472, -5311437055462369280, | |
| /* [ 52] */ -5342038772315650560, -5371315728843297024, -5399331404632512768, -5426144845448965120, | |
| /* [ 56] */ -5451811038519422464, -5476381248265593088, -5499903320558339072, -5522421955752311296, | |
| /* [ 60] */ -5543978956085263616, -5564613449659060480, -5584362093436146432, -5603259257517428736, | |
| /* [ 64] */ -5621337193070986240, -5638626184974132224, -5655154691220933888, -5670949470294763008, | |
| /* [ 68] */ -5686035697601807872, -5700437072199152384, -5714175914219812352, -5727273255295220992, | |
| /* [ 72] */ -5739748920271997440, -5751621603810412032, -5762908939773946112, -5773627565915007744, | |
| /* [ 76] */ -5783793183152377600, -5793420610475628544, -5802523835894661376, -5811116062947570176, | |
| /* [ 80] */ -5819209754516120832, -5826816672854571776, -5833947916825278208, -5840613956570608128, | |
| /* [ 84] */ -5846824665591763456, -5852589350491075328, -5857916778480726528, -5862815203334800384, | |
| /* [ 88] */ -5867292388935742464, -5871355631762284032, -5875011781262890752, -5878267259039093760, | |
| /* [ 92] */ -5881128076579883520, -5883599852028851456, -5885687825288565248, -5887396872144963840, | |
| /* [ 96] */ -5888731517955042304, -5889695949247728384, -5890294025706689792, -5890529289910829568, | |
| /* [100] */ -5890404977675987456, -5889924026487208448, -5889089083913555968, -5887902514965209344, | |
| /* [104] */ -5886366408898372096, -5884482585690639872, -5882252601321090304, -5879677752995027712, | |
| /* [108] */ -5876759083794175232, -5873497386318840832, -5869893206505510144, -5865946846617024256, | |
| /* [112] */ -5861658367354159104, -5857027590486131456, -5852054100063428352, -5846737243971504640, | |
| /* [116] */ -5841076134082373632, -5835069647234580480, -5828716424754549248, -5822014871949021952, | |
| /* [120] */ -5814963157357531648, -5807559211080072192, -5799800723447229952, -5791685142338073344, | |
| /* [124] */ -5783209670985158912, -5774371264582489344, -5765166627072226560, -5755592207057667840, | |
| /* [128] */ -5745644193442049280, -5735318510777133824, -5724610813433666560, -5713516480340333056, | |
| /* [132] */ -5702030608556698112, -5690148005851018752, -5677863184109371904, -5665170350903313408, | |
| /* [136] */ -5652063400924580608, -5638535907000141312, -5624581109999480320, -5610191908627599872, | |
| /* [140] */ -5595360848093632768, -5580080108034218752, -5564341489875549952, -5548136403221394688, | |
| /* [144] */ -5531455851545399296, -5514290416593586944, -5496630242226406656, -5478465016761742848, | |
| /* [148] */ -5459783954986665216, -5440575777891777024, -5420828692432397824, -5400530368638773504, | |
| /* [152] */ -5379667916699401728, -5358227861294116864, -5336196115274292224, -5313557951078385920, | |
| /* [156] */ -5290297970633451520, -5266400072915222272, -5241847420214015744, -5216622401043726592, | |
| /* [160] */ -5190706591719534080, -5164080714589203200, -5136724594099067136, -5108617109269313024, | |
| /* [164] */ -5079736143458214912, -5050058530461741312, -5019559997031891968, -4988215100963582976, | |
| /* [168] */ -4955997165645491968, -4922878208652041728, -4888828866780320000, -4853818314258475776, | |
| /* [172] */ -4817814175855180032, -4780782432601701888, -4742687321746719232, -4703491227581444608, | |
| /* [176] */ -4663154564978699264, -4621635653358766336, -4578890580370785792, -4534873055659683584, | |
| /* [180] */ -4489534251700611840, -4442822631898829568, -4394683764809104128, -4345060121983362560, | |
| /* [184] */ -4293890858708922880, -4241111576153830144, -4186654061692619008, -4130446006804747776, | |
| /* [188] */ -4072410698657718784, -4012466683838401024, -3950527400305017856, -3886500774061896704, | |
| /* [192] */ -3820288777467837184, -3751786943594897664, -3680883832433527808, -3607460442623922176, | |
| /* [196] */ -3531389562483324160, -3452535052891361792, -3370751053395887872, -3285881101633968128, | |
| /* [200] */ -3197757155301365504, -3106198503156485376, -3011010550911937280, -2911983463883580928, | |
| /* [204] */ -2808890647470271744, -2701487041141149952, -2589507199690603520, -2472663129329160192, | |
| /* [208] */ -2350641842139870464, -2223102583770035200, -2089673683684728576, -1949948966090106880, | |
| /* [212] */ -1803483646855993856, -1649789631480328192, -1488330106139747584, -1318513295725618176, | |
| /* [216] */ -1139685236927327232, -951121376596854784, -752016768184775936, -541474585642866432, | |
| /* [220] */ -318492605725778432, -81947227249193216, 169425512612864512, 437052607232193536, | |
| /* [224] */ 722551297568809984, 1027761939299714304, 1354787941622770432, 1706044619203941632, | |
| /* [228] */ 2084319374409574144, 2492846399593711360, 2935400169348532480, 3416413484613111552, | |
| /* [232] */ 3941127949860576256, 4515787798793437952, 5147892401439714304, 5846529325380406016, | |
| /* [236] */ 6622819682216655360, 7490522659874166016, 8466869998277892096, 8216968526387345408, | |
| /* [240] */ 4550693915488934656, 7628019504138977280, 6605080500908005888, 7121156327650272512, | |
| /* [244] */ 2484871780331574272, 7179104797032803328, 7066086283830045440, 1516500120817362944, | |
| /* [248] */ 216305945438803456, 6295963418525324544, 2889316805630113280, -2712587580533804032, | |
| /* [252] */ 6562498853538167040, 7975754821147501312, -9223372036854775808, -9223372036854775808, | |
| ] | |
| var rng: UniformRandomProvider | |
| private var exponential: Exponential | |
| init(rng: UniformRandomProvider) { | |
| self.rng = rng | |
| self.exponential = Exponential(rng: rng) | |
| } | |
| mutating func sample() -> Double { | |
| let xx = rng.nextUniform() | |
| let i = Int(truncatingIfNeeded: xx) & Self.maskInt8 | |
| if i < Self.I_MAX { | |
| return Self.X[i] * Double(Int64(bitPattern: xx)) | |
| } | |
| return edgeSample(initial: xx) | |
| } | |
| private mutating func edgeSample(initial xx: UInt64) -> Double { | |
| let u1 = xx & UInt64(Int64.max) | |
| let signBit = Double(((xx >> 62) & 0x2)) - 1.0 | |
| let j = selectRegion() | |
| var x: Double | |
| if j > Self.J_INFLECTION { | |
| var currentU1 = u1 | |
| while true { | |
| x = Self.interpolate(Self.X, j: j, u: currentU1) | |
| let u = randomInt63() | |
| let uDistance = Int64(truncatingIfNeeded: u) - Int64(truncatingIfNeeded: currentU1) | |
| if uDistance >= 0 { | |
| break | |
| } | |
| if uDistance >= Self.CONVEX_E_MAX && | |
| Self.interpolate(Self.Y, j: j, u: u) < exp(-0.5 * x * x) { | |
| break | |
| } | |
| currentU1 = randomInt63() | |
| } | |
| return signBit * x | |
| } else if j < Self.J_INFLECTION { | |
| if j == 0 { | |
| var value: Double | |
| repeat { | |
| value = Self.ONE_OVER_X_0 * exponential.sample() | |
| } while exponential.sample() < 0.5 * value * value | |
| x = value + Self.X_0 | |
| } else { | |
| var currentU1 = u1 | |
| while true { | |
| let u = randomInt63() | |
| let u2 = currentU1 < u ? u : currentU1 | |
| currentU1 = currentU1 < u ? currentU1 : u | |
| x = Self.interpolate(Self.X, j: j, u: currentU1) | |
| if u2 &- currentU1 > Self.CONCAVE_E_MAX || | |
| Self.interpolate(Self.Y, j: j, u: u2) < exp(-0.5 * x * x) { | |
| break | |
| } | |
| currentU1 = randomInt63() | |
| } | |
| } | |
| return signBit * x | |
| } else { | |
| var currentU1 = u1 | |
| while true { | |
| x = Self.interpolate(Self.X, j: j, u: currentU1) | |
| if Self.interpolate(Self.Y, j: j, u: randomInt63()) < exp(-0.5 * x * x) { | |
| break | |
| } | |
| currentU1 = randomInt63() | |
| } | |
| return signBit * x | |
| } | |
| } | |
| private mutating func selectRegion() -> Int { | |
| let x = rng.nextUniform() | |
| let j = Int(truncatingIfNeeded: x) & Self.maskInt8 | |
| let alias = Int(Self.MAP[j]) & Self.maskInt8 | |
| let signedX = Int64(bitPattern: x) | |
| return signedX >= Self.IPMF[j] ? alias : j | |
| } | |
| } | |
| } |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment