Did anyone implement Ziggurat algorithm already?
I converted one from GSL into C++. In my tests I paired Mersenne Twister with Acklam's approximation, boost's inverse CDF and Ziggurat to generate standard normal distribution. Ziggurat method was 3 times faster than Acklam's approximation and 7 times faster than using boost's inverse CDF. GSL is licensed under GPL v3, but I adopted only "matematical formula", the implementation follows the original Ziggurat paper. My implementation is hardwired to use Mersenne Twister directly and doesn't rely on GSL. GSL version uses smaller look up table and exponential distribution for the tails. I might create my own lookup table, but it would be exactly the same. The math formula is from the published paper, not from GSL. Given that the only thing I took from GSL is precomputed look up table, is my implementation covered by GPL? What do you say? ------------------------------------------------------------------------------ ThinkGeek and WIRED's GeekDad team up for the Ultimate GeekDad Father's Day Giveaway. ONE MASSIVE PRIZE to the lucky parental unit. See the prize list and enter to win: http://p.sf.net/sfu/thinkgeek-promo _______________________________________________ QuantLib-dev mailing list [hidden email] https://lists.sourceforge.net/lists/listinfo/quantlib-dev |
On Tue, Jun 15, 2010 at 5:07 PM, Kakhkhor Abdijalilov
<[hidden email]> wrote: > In my tests I paired Mersenne > Twister with Acklam's approximation, boost's inverse CDF and Ziggurat > to generate standard normal distribution. for those who wonder Acklam's approximation is the (best) inverse CDF implemented in QuantLib > Ziggurat method was 3 times > faster than Acklam's approximation and 7 times faster than using > boost's inverse CDF. ok, so Ziggurat is interesting, even if it is only for pseudo-random and not for quasi-random > GSL version uses smaller look up > table and exponential distribution for the tails. smaller than what ? what is this lookup table ? > I might create my > own lookup table, but it would be exactly the same. of course, if it is the inverse cdf lookup table it's no problem to generate it again (even with Acklam improved by using one iteration of Halley's rational method in order to obtain full machine precision) > GSL is licensed under GPL v3, but I adopted only > "matematical formula", the implementation follows the original > Ziggurat paper. My implementation is hardwired to use Mersenne Twister > directly and doesn't rely on GSL. > [...] The math formula > is from the published paper, not from GSL. Given that the only thing I > took from GSL is precomputed look up table, is my implementation > covered by GPL? I would recompute the lookup table and then it's not covered by GPL ciao -- Nando ------------------------------------------------------------------------------ ThinkGeek and WIRED's GeekDad team up for the Ultimate GeekDad Father's Day Giveaway. ONE MASSIVE PRIZE to the lucky parental unit. See the prize list and enter to win: http://p.sf.net/sfu/thinkgeek-promo _______________________________________________ QuantLib-dev mailing list [hidden email] https://lists.sourceforge.net/lists/listinfo/quantlib-dev |
OK, it is done. I submitted my implementation to SourceForge patch
system (Ziggurat Algorithm by renorm). The underlying RNG is MersenneTwisterUniformRng (MT19937). In my tests Ziggurat was about 2-3 times faster than MT19937+Acklam and about 7 times faster than MT19937+boot's quantile function. I tested Ziggurat with Diehard battery and everything looked OK. I also did some changes to MersenneTwisterUniformRng and InverseCumulativeNormal (Acklam's method). On my platform MersenneTwisterUniformRng is slightly faster than boost's mt19937 and considerable faster that the old MersenneTwisterUniformRng, especially in debug mode. Apparently, using STL vector instead of plain array makes a difference in this case. InverseCumulativeNormal features static std_ic (standard normal inverse cumulative) method. It computes inverse CDF with mean=0 and sigma=0. Using std_ic instead of operator() saves us 2 floating point operation. The performance difference is quite noticeable. Regards, Kakhkhor Abdijalilov. ------------------------------------------------------------------------------ ThinkGeek and WIRED's GeekDad team up for the Ultimate GeekDad Father's Day Giveaway. ONE MASSIVE PRIZE to the lucky parental unit. See the prize list and enter to win: http://p.sf.net/sfu/thinkgeek-promo _______________________________________________ QuantLib-dev mailing list [hidden email] https://lists.sourceforge.net/lists/listinfo/quantlib-dev |
On Wed, 2010-06-16 at 18:22 -0500, Kakhkhor Abdijalilov wrote:
> OK, it is done. I submitted my implementation to SourceForge patch > system Ok, thanks. I'll have a look. Luigi -- A debugged program is one for which you have not yet found the conditions that make it fail. -- Jerry Ogdin ------------------------------------------------------------------------------ ThinkGeek and WIRED's GeekDad team up for the Ultimate GeekDad Father's Day Giveaway. ONE MASSIVE PRIZE to the lucky parental unit. See the prize list and enter to win: http://p.sf.net/sfu/thinkgeek-promo _______________________________________________ QuantLib-dev mailing list [hidden email] https://lists.sourceforge.net/lists/listinfo/quantlib-dev |
I implemented it from the scratch. No copyright problems, whatsoever.
It follows the recipe from Marsaglia and Tsang's paper, except the sampling from the tail. My implementation uses the conditional distribution to sample from the tail of normal density. Regards, Kakhkhor Abdijalilov. ------------------------------------------------------------------------------ Beautiful is writing same markup. Internet Explorer 9 supports standards for HTML5, CSS3, SVG 1.1, ECMAScript5, and DOM L2 & L3. Spend less time writing and rewriting code and more time creating great experiences on the web. Be a part of the beta today. http://p.sf.net/sfu/beautyoftheweb _______________________________________________ QuantLib-dev mailing list [hidden email] https://lists.sourceforge.net/lists/listinfo/quantlib-dev |
Free forum by Nabble | Edit this page |