3
\$\begingroup\$

This is a very specific question mixing up stochastic knowledge and VBA skills. So very exciting!

I'm trying to compare several methods for generating standard, normally distributed numbers given a source of uniformly distributed random numbers. Therefore I'm implementing the Box Muller Algorithm, Ziggurat Algorithm and Ratio of Uniforms Algorithm. Every single implementation works great in terms of generating a clean standard, normally distribution. (checked by Shapiro-Wilk-Test).

What I want to find out: which is the quickest method?

Testing every single program with a total of 10^7 generated numbers these are the run times:

Box-Muller: 3,7 seconds
Ziggurat: 1,28 seconds
Ratio of Uniforms: 10,77 seconds

Actually I am very happy about those readings, because I didn't expect it to be that fast. Of course the run time of every single method also depends on my programming skills and VBA knowledge.

My problem: after doing some research I found out that the Ratio of Uniforms Algorithm should be the quickest (about 3 to 4 times quicker than Box Muller). This information just leans on this stack:

enter image description here

I am curious if this is just a wrong claim of this user or (what I do expect more) if my code is not perfectly implemented. Therefore I'll post my code and hope someone could help me with my question, if my code is just not good enough or if the Ratio of Uniforms just doesn't work that quick as mentioned.

Sub RatioUniforms()

Dim x(10000000) As Double
Dim passing As Long
Dim amount As Long: amount = 10000000
Dim u1 As Double
Dim u2 As Double
Dim v2 As Double

    Do While passing <= amount
        Do
            u1 = Rnd    'rnd= random number(0,1)
        Loop Until u1 <> 0  'u1 musn't become 0
        v2 = Rnd
        u2 = (2 * v2 - 1) * (2 * exp(-1)) ^ (1 / 2)
        If u1 ^ 2 <= exp(-1 / 2 * u2 ^ 2 / u1 ^ 2) Then
            x(passing) = u2 / u1
            passing = passing + 1
        End If
    Loop

 End Sub

Thank you very much helping me on this topic. Maybe some of you have tried those algorithms in VBA or other languages and can help me with there experience about the run time? If you need something else to know about my other implementations just let me know. Have a great day!

\$\endgroup\$
3
\$\begingroup\$

There are a few optimizations that can be made to the code to speed it up, mainly to do with how some mathematical operations are performed in VBA. Not sure if I have it implemented 100% accurately so please review for accuracy. Also, I'm sure there are more optimizations to be had, but this is hopefully a starting point for further conversation.

A lot of my changes may be PC specific. It may ultimately depend on your CPU and instruction set available. On my computer, this runs in about 2.5 seconds.

A list of changes:

  1. I replaced all instances of raising something to a power of two, instead, I just multiplied the item by itself.
  2. I pre-computed this part Sqr((2 * Exp(-1))) as it appears to always be the same, so it isn't calculated for each loop and put it into a constant.
  3. I removed the variable v2, it wasn't really needed, you just needed to introduce Rnd in one more spot
  4. General cleanup of the code, and renamed a few variables for clarity

Code

Sub RatioUniforms()
    Const NumberOfIterations        As Long = 10000000
    Const u2CalculationSecondHalf   As Double = 0.857763884960707 'Caching this part Sqr((2 * Exp(-1)))
    Dim Results(NumberOfIterations) As Double
    Dim PassCounter                 As Long
    Dim u1                          As Double 'Define a better name if possible
    Dim u2                          As Double 'Define a better name if possible
    Dim MyTimer                     As Double

    MyTimer = Timer
    Do While PassCounter <= NumberOfIterations
        Do: u1 = Rnd: Loop Until u1 > 0
        u2 = (2 * Rnd - 1) * u2CalculationSecondHalf

        If u1 * u1 <= Exp(-1 / 2 * (u2 * u2) / (u1 * u1)) Then
            Results(PassCounter) = u2 / u1
            PassCounter = PassCounter + 1
        End If

    Loop
    Debug.Print "Process took: " & Timer - MyTimer

 End Sub
\$\endgroup\$
  • \$\begingroup\$ Thank you very much for your adviece! Now it runs also 3,5 seconds. Stil comparing to Box Muller it is not 3x faster... Do you know about this topic? Is Ratio of Uniforms really that much faster? \$\endgroup\$ – J.schmidt Jan 9 at 14:08
  • \$\begingroup\$ @J.schmidt. Probably theoretically, but it very much depends on the implementation. I've updated my answer slightly, got it improved to about 2.5 seconds. If you are looking for pure speed, VBA probably not the right tool for the job. \$\endgroup\$ – Ryan Wildry Jan 9 at 14:11
  • \$\begingroup\$ What else did you change to improve it to 2,5 seconds? I'm aware that in terms of run time VBA is not the best too. Stil, my task is to see which algorithm is the fastest using VBA. \$\endgroup\$ – J.schmidt Jan 9 at 15:08
  • 2
    \$\begingroup\$ I cached this calculation Sqr((2 * Exp(-1))) as it appears to be deterministic and put that into a Constant Variable. \$\endgroup\$ – Ryan Wildry Jan 9 at 15:09
  • \$\begingroup\$ The only idea I had, was to simplify the test u1 * u1 <= Exp(-1 / 2 * (u2 * u2) / (u1 * u1)). It doesn't match ~30% of time, so there are quite a bit of wasted cycles here. Not comfortable enough to know how to simplify this test. The one thing that stuck out, is u1*u1 is on both sides. If there was a way to simplify this expression, it might run even faster. I thought something like Exp(-1 / 2 * (u2 * u2) / (u1 * u1)) / u1 * u1 <= 1 would work, but not qualified to say for sure \$\endgroup\$ – Ryan Wildry Jan 9 at 16:44

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.