I've tried to find documentation on how the Delphi V7 random number
generator is implemented but so far failed. It does fail most DIEHARD
tests badly but it would still be nice to know what kind of generator it
is. I suspect some LCG but is too short on time to try to figure out if
that is the case and if so, what the parameters are.
Anyone knows?
IMO, it's a pseudo randomizer. It generates a random number based on a
base number, the random seed. The Random function does the calculation
and it'll set the random seed value based on the previous random seed.
The Randomize function on the other hand, set the random seed from the
computer timer counter or the CPU clock counter in later version.
Analyzing both function will get you the idea of how the random number
is generated.
Yes, but what family of RNGs does it belong to and in case it is some
common family (such as linear congruential), what are the parameters? It
fails pretty much every DIEHARD test for pretty much all bits, both hi
and low, even if one only uses truncated outputs from the PRNG.
The solution was for us to make an implementation of the KISS PRNG in
Delphi so we don't have any problems generating nice looking
pseudo-randomness. We would be helped by knowing the exact
implementation of the original PRNG since it could explain some odd
behaviours from the original program.
I hope converting a program from one language to another is just a small
step of your overall project. By itself, that's not a very interesting
master's task.
> The original program does a simulation
> relying heavily on Random(arg) and Random.
>
> I've tried to find documentation on how the Delphi V7 random number
> generator is implemented but so far failed. It does fail most DIEHARD
> tests badly but it would still be nice to know what kind of generator it
> is. I suspect some LCG but is too short on time to try to figure out if
> that is the case and if so, what the parameters are.
You have Delphi, so you can look at the source code. It's in System.pas.
The implementation is subject to change between Delphi releases. I'll
describe the Delphi 2005 implementation:
Take RandSeed, multiply by 0x08088405, and truncate to 32 bits.
Increment the result and store that as RandSeed. Multiply the value with
the input range parameter for Random, and return the upper 32 bits.
The floating-point version is the same, but instead of multiplying the
new RandSeed with the integer range value, multiply by 2^-32.
--
Rob
The conversion itself is a minor part indeed. The original program's
functionality is a subset of ours, but we need to, as closely as
possible, emulate the behaviour of the old program.
[snip]
> You have Delphi, so you can look at the source code. It's in System.pas.
> The implementation is subject to change between Delphi releases. I'll
> describe the Delphi 2005 implementation:
>
> Take RandSeed, multiply by 0x08088405, and truncate to 32 bits.
> Increment the result and store that as RandSeed. Multiply the value with
> the input range parameter for Random, and return the upper 32 bits.
>
> The floating-point version is the same, but instead of multiplying the
> new RandSeed with the integer range value, multiply by 2^-32.
Thanks a lot!
I twiced analyzed and tried to reproduce it in simply pascal/delphi code.
The first time unsuccessfull.
The second time successfull.
I even made a 64 bit version if I recall correctly ;)
Google Skybuck, Random Number Generator
Let me know if you find it... if not I might repost it :):):)
Bye,
Skybuck ;)
>> The original program does a simulation relying heavily on
>>Random(arg) and Random.
>> I've tried to find documentation on how the Delphi V7 random number
>>generator is implemented but so far failed. It does fail most DIEHARD
>>tests badly but it would still be nice to know what kind of generator
>>it is. I suspect some LCG but is too short on time to try to figure
>>out if that is the case and if so, what the parameters are.
>
>You have Delphi, so you can look at the source code. It's in
>System.pas. The implementation is subject to change between Delphi
>releases. I'll describe the Delphi 2005 implementation:
>
>Take RandSeed, multiply by 0x08088405, and truncate to 32 bits.
>Increment the result and store that as RandSeed. Multiply the value
>with the input range parameter for Random, and return the upper 32
>bits.
>
>The floating-point version is the same, but instead of multiplying the
>new RandSeed with the integer range value, multiply by 2^-32.
Some of that, and consequences, is in my
<URL:http://www.merlyn.demon.co.uk/pas-rand.htm>.
<URL:http://www.merlyn.demon.co.uk/programs/mathutys.pas> has a 64-bit
version, slightly tested in
<URL:http://www.merlyn.demon.co.uk/programs/mathunit.zip>.
--
(c) John Stockton, nr London, UK. ?@merlyn.demon.co.uk Turnpike v6.05 MIME.
Web - FAQqish topics, acronyms & links;
Astro stuff via astron-1.htm, gravity0.htm ; quotings.htm, pascal.htm, etc.
No Encoding. Quotes before replies. Snip well. Write clearly. Don't Mail News.
[snip]
> Some of that, and consequences, is in my
> <URL:http://www.merlyn.demon.co.uk/pas-rand.htm>.
> <URL:http://www.merlyn.demon.co.uk/programs/mathutys.pas> has a 64-bit
> version, slightly tested in
> <URL:http://www.merlyn.demon.co.uk/programs/mathunit.zip>.
Thanks. From the assembly version I translated it to this (C), matching
the output from Delphi;
---
static unsigned int randseed;
int rnd(int m) {
randseed = randseed * 0x08088405ULL + 1;
return (int) ((randseed * (unsigned long long) m) >> 32);
}
---
We did end up using KISS for our Delphi needs, replacing the original
PRNG to make sure it didn't cause any strange behaviour. Our KISS
implementation looks like this;
---
unit kiss;
interface
type
RNGState = class(TObject)
x, y, z, w : LongWord;
end;
MyRandom = class(TObject)
public
Constructor Init(Seed : LongWord);
function GetInt(N : Integer) : Integer;
function GetFloat : Extended;
function GetLW() : LongWord;
private
RState : RNGState;
end;
implementation
{$Q-}
function KISSRandom(RState : RNGState) : LongWord;
begin
RState.x := RState.x * 69069 + 1234567;
RState.y := RState.y xor RState.y shl 13;
RState.y := RState.y xor RState.y shr 17;
RState.y := RState.y xor RState.y shl 5;
RState.z := 65184 * (RState.z and 65535) + (RState.z shr 16);
RState.w := 64860 * (RState.w and 65535) + (RState.w shr 16);
KISSRandom := (((RState.z shl 16) + RState.w) xor RState.x) + RState.y;
end;
procedure SeedKISS(RState : RNGState; Seed : LongWord);
var
I : Integer;
begin
// x has period 2^32 so it can be seeded with anything and still have
// full period.
RState.x := Seed;
// y has period 2^32 - 1, having cycles of length {{1}, {2^32-1}}
if Seed = 0 then
RState.y := 1
else
RState.y := Seed;
// For the MWC generators, the carry must not be equal to the
// multipliers - 1.
RState.z := (1 shl 16) + (Seed shr 16);
RState.w := (1 shl 16) + (Seed And $ffff);
// Now run the generator for a few steps;
for I := 0 to 31 do begin
KISSRandom(RState);
end;
// Done with seeding.
end;
Constructor MyRandom.Init(Seed : LongWord);
begin
RState := RNGState.Create();
SeedKISS(RState, Seed);
end;
function MyRandom.GetInt(N : Integer) : Integer;
var
Bits, Val : LongWord;
begin
// Taken from the java.util.Random manpage, which in turn borrows from
// Warren's "Hacker's Delight".
{
if (N and -N) = N then begin // i.e., n is a power of 2
GetInt := LongWord((UInt64(N) * KISSRandom(RState)) shr 32);
Exit;
end;
}
repeat
Bits := KISSRandom(RState);
Val := Bits mod N;
until Bits - Val + (N - 1) >= 0;
GetInt := Val;
end;
function MyRandom.GetFloat : Extended;
begin
GetFloat := KISSRandom(RState) / UInt64($100000000);
end;
function MyRandom.GetLW : LongWord;
begin
GetLW := KISSRandom(RState);
end;
end.
---
We have very little experience with Pascal/Delphi which may be very
obvious. :-)
Anyway, GetInt(N) should generate uniform numbers on [0, N) and GetFloat
(almost) uniform floats on [0, 1). The generator has a period close to
2^126 (2^15 * 64860 - 1) * (2^15 * 65184 - 1) * 2^32 * (2^32 - 1). See
Marsaglia's 1999 Usenet posting;
http://oldmill.uchicago.edu/~wilder/Code/random/Papers/Marsaglia_1999.html
and
http://www.evensen.org/marsaglia/mwc1.ps (Not sure where I found that
paper, it seems to not be available anyplace else public anymore.)
> if (N and -N) = N then begin // i.e., n is a power of 2
I think N and (N-1) = 0 is equivalent; it might be
(insignificantly) quicker.
There is a small problem available when converting a 64-bit random
unsigned integer to an IEEE Double X such that 0.0 <= X < 1.0 . If one
just float-divides by 2^64, the absolute density of results increases
for each of the first about 9 times that the result exponent decreases -
that's because the float has a 53-bit mantissa. Conversion to 64-bit
mantissa extended should be OK, but the problem will become included if
that extended is rounded to Double. The answer seems to be to mask the
unsigned integer down to 53 bits before division - but think about which
bits to mask out for best results.
--
(c) John Stockton, nr London UK. ?@merlyn.demon.co.uk BP7, Delphi 3 & 2006.
<URL:http://www.merlyn.demon.co.uk/> TP/BP/Delphi/&c., FAQqy topics & links;
<URL:http://www.bancoems.com/CompLangPascalDelphiMisc-MiniFAQ.htm> clpdmFAQ;
NOT <URL:http://support.codegear.com/newsgroups/>: news:borland.* Guidelines