While numerically testing a solution for the
February 2007 IBM Ponder This challenge
I had to generate a large n×2 Array of random floats
uniformly distributed in [0,1).

My first attempt was to use *structured flavors* in the
`RandomTools`

package:

picks[4] := proc(n::posint)
uses RandomTools;
Array(Generate(list(list(float(range=0..1,'method=uniform'),2),n))
, datatype=float[8]);
end proc:

While this worked, it was rather slow. I then noticed that
there is a `listlist`

flavor in the
`RandomTools`

package. Using it gave a significant
improvement (about 2×) over using nested `list`

flavors.

picks[3] := proc(n::posint)
uses RandomTools;
Array(Generate(listlist(float(range=0..1,'method=uniform'),n,2))
, 'datatype=float[8]'
);
end proc:

By specifying a 2×n `listlist`

structure and
using the `transpose`

option to the
`Array`

constructor this could be slightly improved:

picks[2] := proc(n::posint)
uses RandomTools;
Array(Generate(listlist(float(range=0..1,'method=uniform'),2,n))
, 'datatype=float[8]'
, 'transpose=true'
);
end proc:

While the structured flavors are convenient, they are not as
fast as I'd like. My next attempt was to use the
`GenerateFloat64`

procedure in the
`MersenneTwister`

subpackage of
`RandomTools`

. It calls (I believe) compiled code
that is significantly faster than the random number generator
used to generate structured flavors. Unlike the structured flavor
approach, `GenerateFloat64`

returns only a singly random float
that is (thankfully) uniformly distributed on [0,1) so we'll have to
use an appropriate constructor to initialize the Array.

picks[1.3] := proc(n::posint)
uses RandomTools;
Array(1..n, 1..2, (i,j) -> MersenneTwister:-GenerateFloat64()
, 'datatype=float[8]'
);
end proc:

That is significantly faster (about 5× with `n=50000`

)
than the previous method. As you may have guessed by the numbering
scheme, there are still a few improvements.
Each is somewhat faster than the next.
picks[1.2] := proc(n::posint)
uses RandomTools;
local i;
Array([seq([MersenneTwister:-GenerateFloat64()
, MersenneTwister:-GenerateFloat64()]
,i=1..n)
]
, 'datatype=float[8]'
);
end proc:
picks[1.1] := proc(n::posint)
uses RandomTools;
local i;
Array([NULL
, [seq(MersenneTwister:-GenerateFloat64(), i=1..n)]
, [seq(MersenneTwister:-GenerateFloat64(), i=1..n)]
]
, 'transpose = true'
, 'datatype=float[8]'
);
end proc:

Jacques Carette suggests the following trick, which directly calls the compiled procedure
used by GenerateFloat64.

picks[1.0] := proc(n::posint) local t;
uses RandomTools; kernelopts(opaquemodules=false):
t := (a,b) -> MersenneTwister:-MTKernelInterface(3);
Array(1..n, 1..2
, t
, 'datatype'=float[8]
);
end proc:

Here is a minor improvement. First, I restored the setting of `opaquemodules`

. I also used the faster technique of passing a `listlist`

structure rather than using functional initializer. Finally, to generate the lists quickly, two `seq`

commands were used and the structure transposed, using the transpose option to Array.

picks[0.8] := proc(n::posint)
local i,rnd,opacity;
uses RandomTools;
opacity := kernelopts('opaquemodules'=false);
rnd := MersenneTwister:-MTKernelInterface;
kernelopts('opaquemodules'=opacity);
Array([NULL
, [seq(rnd(3), i=1..n)]
, [seq(rnd(3), i=1..n)]
]
, 'transpose=true'
, 'datatype=float[8]'
);
end proc:

This, too, can be further improved. For this case it is faster to create an empty Array and then use the `map`

function to initialize random data.

picks[0.7] := proc(n::posint)
local opacity,rnd;
opacity := kernelopts('opaquemodules'=false);
rnd := RandomTools:-MersenneTwister:-MTKernelInterface;
kernelopts('opaquemodules'=opacity);
map(()->rnd(3), Array(1..n,1..2,'datatype=float[8]'));
end proc:

A bit more speed can be gained by using a hack. A number, including a hardware float, that is applied to any argument evaluates to itself. For example, `3.2(0.0) = 3.2`

. We can use this trick to avoid having to evaluate the operator `()->rnd(3)`

, instead we use the unevaluated function `'rnd(3)'`

which, when evaluated, returns a hardware float that is then applied to 0.0 (the initial value of the Array) and so evaluates to itself:

picks[0.65] := proc(n::posint)
local opacity,rnd;
opacity := kernelopts('opaquemodules'=false);
rnd := RandomTools:-MersenneTwister:-MTKernelInterface;
kernelopts('opaquemodules'=opacity);
map('rnd(3)', Array(1..n,1..2,'datatype=float[8]'));
end proc:

Acer points out that `rtable`

has a built-in `frandom`

initializer. While the documentation for this feature is deficient (it doesn't mention whether the generated random numbers are uniformly distributed within the range, or other details) it is orders of magnitude faster than the previous approach.

picks[0.1] := proc(n::posint)
rtable(1..n, 1..2
, frandom(0..1, 1)
, 'subtype = Array'
, 'datatype=float[8]'
);
end proc:

Here is a performance comparison of all the techniques.
and scaling the result.
The fastest procedure is abpit 3 orders of magnitude faster than my first attempt.
All tests were run with Maple 11 on a Linux box.

GGenerate a 50000 x 2 Array of random floats uniformly distributed in [0,1)
----- timing ------- ---------- memory ------------
total (s) maple (s) used (MB) alloc (MB) gctimes proc
--------- --------- --------- ---------- ------- ----
0.02 0.02 0.76 0.81 0 picks[.1]
0.85 0.82 18.32 10.87 1 picks[.65]
1.03 0.92 18.96 10.69 1 picks[.7]
1.14 1.11 20.82 13.19 1 picks[.8]
1.25 1.20 21.97 14.62 1 picks[.85]
1.37 1.37 21.64 11.75 1 picks[.9]
1.50 1.46 21.74 11.75 1 picks[1.0]
1.41 1.37 23.43 14.31 1 picks[1.1]
1.53 1.50 24.12 16.06 1 picks[1.2]
1.68 1.61 23.41 13.25 1 picks[1.3]
7.79 7.58 125.29 39.74 4 picks[2]
8.51 8.36 129.21 41.49 4 picks[3]
16.63 16.01 199.61 43.93 6 picks[4]