BBC BASIC implementation of RANDU

[p_m21987]

Hello,
I've recently been reading about the history of pseudo-random number generators, and I've found RANDU to be particularly fascinating.

So tonight for fun, I implemented RANDU in BBC BASIC and made this little graphics demo with it:

Code: Select all

      randu_seed%%=1

      n%=1600
      DIM p(n%,2), c(n%)
      FOR i%=0 TO n%
        FOR j%=0 TO 2
          p(i%,j%) = (FNrandu-0.5)*700
        NEXT
        c(i%)=1+FNrandu*14
      NEXT

      DIM r(2,2), r2(2,2), a(n%,2)

      MODE 8
      OFF
      ORIGIN 640,512
      *refresh off
      REPEAT
        MOUSE X%,Y%,B%
        a=X%/640*PI
        a2=Y%/-512*PI
        r(0,0)=COSa: r(2,0)=SINa: r(1,1)=1.0: r(0,2)=-SINa: r(2,2)=COSa
        r2(0,0)=1.0: r2(1,1)=COSa2: r2(2,1)=-SINa2: r2(1,2)=SINa2: r2(2,2)=COSa2
        a() = p() . r()
        a() = a() . r2()
        CLG
        FOR i%=0 TO n%
          GCOL 0, c(i%)
          PLOT 69, a(i%,0), a(i%,1)
        NEXT
        *refresh
        WAIT 1
      UNTIL FALSE

      END

      DEF FNrandu
      randu_seed%% = (65539 * randu_seed%%) MOD 2147483648
      = randu_seed%% / 2147483648

EDIT: Added Richard's improvements, thank you!

The demo generates a number of points in 3D space, and lets you rotate them around by moving the mouse.
From some angles it looks like the points are dispersed randomly... but then when you move the perspective a bit, you'll find that the points all lie on 15 planes. Fascinating!

I've really enjoyed looking into this fascinating and amazing function and learning about the history of it tonight. I hope you will too.
Here is the wikipedia which is a very good read: https://en.wikipedia.org/wiki/RANDU

Regards,
PM

[Hated Moron]

I've recently been reading about the history of pseudo-random number generators, and I've found RANDU to be particularly fascinating.

Fascinating indeed, and a cautionary tale. Just in case anybody skimps the background reading, it's perhaps worth emphasising that "RANDU is widely considered to be one of the most ill-conceived random number generators ever designed"!

This is what it says in the article Notes on the use of RND at the BBC BASIC Wiki: "One thing not to attempt is to write your own algorithm, unless you're a top-notch mathematician!".

So tonight for fun, I implemented RANDU in BBC BASIC and made this little graphics demo with it:

Nice visualisation, and an excellent use of the dot-product operator too!

To avoid using excessive CPU time if run in BBC BASIC for Windows I'd advise changing the WAIT to WAIT 1.

I'd also suggest changing MODE 3 to MODE 8, which will allow it to run nicely in Matrix Brandy as well.

[p_m21987]

Fascinating indeed, and a cautionary tale. Just in case anybody skimps the background reading, it's perhaps worth emphasising that "RANDU is widely considered to be one of the most ill-conceived random number generators ever designed"!

Yes! It's amazing and fascinating that this function was once so widely deployed.

Nice visualisation, and an excellent use of the dot-product operator too!

To avoid using excessive CPU time if run in BBC BASIC for Windows I'd advise changing the WAIT to WAIT 1.

I'd also suggest changing MODE 3 to MODE 8, which will allow it to run nicely in Matrix Brandy as well.

Ahh, thank you! I appreciate it a lot.
I've edited the program listing to include these improvements now too.

Regards, PM

[Hated Moron]

Yes! It's amazing and fascinating that this function was once so widely deployed.

Can I assume that BBC BASIC's built-in RND(1) function doesn't exhibit any obvious weaknesses when subjected to the 3D test? I note that RANDU's sequence length is 2^29 (536870912) whereas RND's sequence length is 2^33-1 (8589934591). Neither of those is long enough for critical applications, such as shuffling a pack of cards.

[p_m21987]

Yes! It's amazing and fascinating that this function was once so widely deployed.

Can I assume that BBC BASIC's built-in RND(1) function doesn't exhibit any obvious weaknesses when subjected to the 3D test? I note that RANDU's sequence length is 2^29 (536870912) whereas RND's sequence length is 2^33-1 (8589934591). Neither of those is long enough for critical applications, such as shuffling a pack of cards.

Testing it with 19200 x,y,z points generated with RND(1), it looked perfectly random to my eyes from every possible angle I could rotate the cube into with the mouse.
Shuffling a pack of cards is an interesting everyday life example. The number of permutations is quite astounding when you think about it.

[Hated Moron]

Shuffling a pack of cards is an interesting everyday life example. The number of permutations is quite astounding when you think about it.

The Wiki article says it's around 8*10^67. A couple of PRNGs that do have a long enough sequence length are listed here although as it says a much better method of seeding them would be required.

[Hated Moron]

Testing it with 19200 x,y,z points generated with RND(1), it looked perfectly random to my eyes from every possible angle I could rotate the cube into with the mouse.

I note that the relevant Wikipedia article says that the 'spectral test' can only be applied to Linear Congruential Generators (LCG), so as it turns out I was wasting your time suggesting that it be used on BBC BASIC's RND(1), which is a Linear Feedback Shift Register (LFSR). Sorry.

I have every reason to believe that the LFSR used by BBC BASIC (right back to the 1981 6502 version) has excellent properties as a PRNG, with its only major shortcoming being the sequence length (2^33-1) which although more than some others limits the applications to which it can be put.

There's a long thread at Stardot discussing almost every conceivable aspect of BBC BASIC's random number generation. ' Deleted User 9295' is me, before I was expelled.

[Hated Moron]

Here's a version which should run in ARM BBC BASIC V (64-bit integer variables eliminated):

Code: Select all

   10 randu_seed%=1
   20
   30 n%=1600
   40 DIM p(n%,2), c(n%)
   50 FOR i%=0 TO n%
   60   FOR j%=0 TO 2
   70     p(i%,j%) = (FNrandu-0.5)*700
   80   NEXT
   90   c(i%)=1+FNrandu*14
  100 NEXT
  110
  120 DIM r(2,2), r2(2,2), a(n%,2)
  130
  140 MODE 8
  150 OFF
  160 ORIGIN 640,512
  170 REPEAT
  180   MOUSE X%,Y%,B%
  190   a=X%/640*PI
  200   a2=Y%/-512*PI
  210   r(0,0)=COSa: r(2,0)=SINa: r(1,1)=1.0: r(0,2)=-SINa: r(2,2)=COSa
  220   r2(0,0)=1.0: r2(1,1)=COSa2: r2(2,1)=-SINa2: r2(1,2)=SINa2: r2(2,2)=COSa2
  230   a() = p() . r()
  240   a() = a() . r2()
  250   CLS
  260   FOR i%=0 TO n%
  270     GCOL 0, c(i%)
  280     PLOT 69, a(i%,0), a(i%,1)
  290   NEXT
  300   WAIT
  310 UNTIL FALSE
  320
  330 END
  340
  350 DEF FNrandu
  360 LOCAL temp%
  370 temp% = (randu_seed% << 16) AND &7FFFFFFF
  380 temp% = (temp% + randu_seed%) AND &7FFFFFFF
  390 temp% = (temp% + randu_seed%) AND &7FFFFFFF
  400 randu_seed% = (temp% + randu_seed%) AND &7FFFFFFF
  410 = randu_seed% / 2147483648

[p_m21987]

Testing it with 19200 x,y,z points generated with RND(1), it looked perfectly random to my eyes from every possible angle I could rotate the cube into with the mouse.

I note that the relevant Wikipedia article says that the 'spectral test' can only be applied to Linear Congruential Generators (LCG), so as it turns out I was wasting your time suggesting that it be used on BBC BASIC's RND(1), which is a Linear Feedback Shift Register (LFSR). Sorry.

No need to be sorry! I think it was still worth trying out regardless for the experience and for the fun of the visualisation.

There's a long thread at Stardot discussing almost every conceivable aspect of BBC BASIC's random number generation. ' Deleted User 9295' is me, before I was expelled.

Very interesting thread, thank you for linking it. I've been reading it. Through reading the thread, I found the C code for the implementation of RND in your BBCSDL, very cool to see and read.
I'm sad you were pushed away from the Stardot forums, I hope you'll be welcomed back. I'm sure the people there would love to see you on the forums again, you're the legendary inventor of BBC BASIC after all.

I have every reason to believe that the LFSR used by BBC BASIC (right back to the 1981 6502 version) has excellent properties as a PRNG, with its only major shortcoming being the sequence length (2^33-1) which although more than some others limits the applications to which it can be put.

I'm really impressed by it and I'm interested in knowing more about it.
After following the stardot link and finding the C code for RND in BBCSDL, I thought it would be neat to be able to use BBC BASIC's RND function in my own C programs.
I hope this is okay with you Richard, let me know if it isn't and I'll delete it if so - I've tried to make a general purpose C adaptation of RND from BBCSDL for my own private use in my own C programming projects.

[Hated Moron]

I'm sure the people there would love to see you on the forums again

Members, maybe. The owner and principal admin of the StarDot forum, definitely not sadly.

I feel bad about 'cheating' my way back onto the Discussion Group - I never post there because it was made very clear that I am not welcome, but they nevertheless see the posts I make here because of the automatic forwarding.

you're the legendary inventor of BBC BASIC after all.

Eek, that wasn't me it was Sophie Wilson. Don't make that mistake at StarDot or the RISC OS forums.

hope this is okay with you Richard, let me know if it isn't and I'll delete it if so

Nothing to do with me. I don't know if Acorn first used that algorithm in BBC BASIC or whether it was previously used in Atom BASIC. Either way it's a standard maximal-length LFSR, as for example can be found listed in this Xilinx application note.