[FFmpeg-devel] [PATCH] lavu/lfg: switch to Ziggurat algorithm for normal random number generation

[Ganesh Ajjanagadde]

On Thu, Mar 10, 2016 at 11:30 PM, James Almer <jamrial at gmail.com> wrote:
> On 3/11/2016 1:16 AM, Ganesh Ajjanagadde wrote:
>>  void av_bmg_get(AVLFG *lfg, double out[2])
>>  {
>> -    double x1, x2, w;
>> -
>> -    do {
>> -        x1 = 2.0 / UINT_MAX * av_lfg_get(lfg) - 1.0;
>> -        x2 = 2.0 / UINT_MAX * av_lfg_get(lfg) - 1.0;
>> -        w  = x1 * x1 + x2 * x2;
>> -    } while (w >= 1.0);
>> -
>> -    w = sqrt((-2.0 * log(w)) / w);
>> -    out[0] = x1 * w;
>> -    out[1] = x2 * w;
>> +    out[0] = ziggurat(lfg);
>> +    out[1] = ziggurat(lfg);
>>  }
>>
>>  #ifdef TEST
>> diff --git a/libavutil/lfg.h b/libavutil/lfg.h
>> index ec90562..6241359 100644
>> --- a/libavutil/lfg.h
>> +++ b/libavutil/lfg.h
>> @@ -52,7 +52,7 @@ static inline unsigned int av_mlfg_get(AVLFG *c){
>>  }
>>
>>  /**
>> - * Get the next two numbers generated by a Box-Muller Gaussian
>> + * Get two numbers generated independently from a standard Gaussian distribution
>
> I don't think a function called after Box-Muller Gaussian should suddenly
> start using a different algorithm to generate random numbers.
>
> How about adding a new av_ziggurat_get() or similarly named function for
> this, and then replace av_bmg_get() usage where convenient?

Fair enough. Personally, I have no idea why a public API function
should name itself after an implementation detail, like Box-Muller,
Ziggurat, ratio of uniforms, inverse cdf, etc. Furthermore, the API is
brittle, the 2 samples again reflects an implementation detail.

For example, not suggesting that we do this right now, but for future
proofing it may be useful to fill up arrays of arbitrary length rather
than just 2. This could help in some SIMD'ing and efficient
utilization of samples. For example, one could use a fast method for
generating a random block, extract from it on a as-needed basis.
Again, gains can be had from vectorization.

BTW, above is not just theoretical but is actually done by one recent
implementation that seems to be > 2x faster than this one:
http://arxiv.org/pdf/1403.6870.pdf.

As such, at the minimum, I suggest naming it as av_gaussian_get or
something like that. As for the signature, no idea but (AVLFG *lfg,
double *out, int len) sounds better, but I am not a good person for
this.

>
>>   * generator using the random numbers issued by lfg.
>>   *
>>   * @param out array where the two generated numbers are placed
>
> _______________________________________________
> ffmpeg-devel mailing list
> ffmpeg-devel at ffmpeg.org
> http://ffmpeg.org/mailman/listinfo/ffmpeg-devel