Bernard Widynski's Middle Square Weyl Sequence RNG (MsWs)

[MrSwiss]

Huston, deltarho[] has a problem!

After giving you time to ponder my last post and, not receiving a reply.
A rhetorical question first: "have you figured it out, by now?"
(obvious answer, seems to be: NO)

And, as proof of concept, my figures from a -5.0 to +5.0 range, with Doubles,
of course and, the value the result referres to (no useless frills attached):

Code: Select all

-5                       4999118
-4.9                     4999118
-4.800000000000001       4999118
-4.700000000000001       4999118
-4.600000000000001       4999118
-4.500000000000002       4999118
-4.400000000000002       10002651
-4.300000000000003       10002651
-4.200000000000003       10002651
-4.100000000000003       10002651
-4.000000000000004       10002651
-3.900000000000004       10002651
-3.800000000000003       10002651
-3.700000000000003       10002651
-3.600000000000003       10002651
-3.500000000000003       10002651
-3.400000000000003       9999430
-3.300000000000003       9999430
-3.200000000000003       9999430
-3.100000000000003       9999430
-3.000000000000003       9999430
-2.900000000000003       9999430
-2.800000000000003       9999430
-2.700000000000002       9999430
-2.600000000000002       9999430
-2.500000000000002       9999430
-2.400000000000002       9997659
-2.300000000000002       9997659
-2.200000000000002       9997659
-2.100000000000002       9997659
-2.000000000000002       9997659
-1.900000000000002       9997659
-1.800000000000002       9997659
-1.700000000000002       9997659
-1.600000000000001       9997659
-1.500000000000001       9997659
-1.400000000000001       10000457
-1.300000000000001       10000457
-1.200000000000001       10000457
-1.100000000000001       10000457
-1.000000000000001       10000457
-0.9000000000000009      10000457
-0.8000000000000009      10000457
-0.700000000000001       10000457
-0.600000000000001       10000457
-0.500000000000001       10000457
-0.400000000000001       10005042
-0.300000000000001       10005042
-0.200000000000001       10005042
-0.100000000000001       10005042
-1.02695629777827e-015   10005042
 0.09999999999999898     10005042
 0.199999999999999       10005042
 0.299999999999999       10005042
 0.399999999999999       10005042
 0.499999999999999       10005042
 0.599999999999999       10000432
 0.699999999999999       10000432
 0.7999999999999989      10000432
 0.8999999999999989      10000432
 0.9999999999999989      10000432
 1.099999999999999       10000432
 1.199999999999999       10000432
 1.299999999999999       10000432
 1.399999999999999       10000432
 1.499999999999999       10000432
 1.599999999999999       10004424
 1.7                     10004424
 1.8                     10004424
 1.9                     10004424
 2                       10004424
 2.1                     10004424
 2.2                     10004424
 2.3                     10004424
 2.4                     10004424
 2.5                     10004424
 2.6                     9997372
 2.7                     9997372
 2.8                     9997372
 2.9                     9997372
 3                       9997372
 3.100000000000001       9997372
 3.200000000000001       9997372
 3.300000000000001       9997372
 3.400000000000001       9997372
 3.500000000000001       9994098
 3.600000000000001       9994098
 3.700000000000001       9994098
 3.800000000000001       9994098
 3.900000000000001       9994098
 4.000000000000001       9994098
 4.100000000000001       9994098
 4.2                     9994098
 4.3                     9994098
 4.4                     9994098
 4.499999999999999       9994098
 4.599999999999999       4999317
 4.699999999999998       4999317
 4.799999999999998       4999317
 4.899999999999998       4999317
 4.999999999999997       4999317

to high: 0   to low: 0   count: 100000000
out of range: 0

press a key to EXIT ...

They show a clearly visible pattern, as oposed to 'integer' data-types.
It also proofes, the Banker's rounding problem, affecting the distribution.

Cheers.

[deltarho[1859]]

@MrSwiss

Your change of tack threw me; or should I say your change of attack. <smile>

You have lost every point on the LongInt discussion and decided to 'have a go' at the Doubles.

I mentioned in another thread sometime ago the problem with FB's generators not providing an integer output forcing us to use RND to get back into the integer domain. Integer => Floating point (via RND) => Integer has issues and it is better to use Integer => Integer.

For decimal ranges FB Help suggests

Code: Select all

Function = Rnd * (last - first) + first

If we look at my opening post we have

Code: Select all

Function = This.Seed/2^64*( Last - First ) + First

This does not have distribution issues.

[deltarho[1859]]

This does not have distribution issues.

Actually, it does.

From FB's manual we have

Code: Select all

Function rnd_range (first As Double, last As Double) As Double
  Function = Rnd * (last - first) + first
End Function

I have never questioned that. In fact, I don't think that anybody has. It has a distribution issue.

However, if we use dafhi's method and then correct that by replacing '+1' with my '-(first>=0)' we get

Code: Select all

Function rnd_range (first As Double, last As Double) As Double
  Function = Rnd*(last - first -(first>=0)) - .5 + first
End Function 

That does not seem to have any distribution issues.

@MrSwiss

Perhaps you, and others can have a look at this because a lot of folks will be using FB's rnd_range.

[MrSwiss]

Your change of tack threw me; or should I say your change of attack.

Now, you are in denial.
I've clearly stated 'the case', see: Function with Doubles

---

I'll look into your 'hack', which needs to be made better readable, anyhow.
Since, it is 'bad practice' IMO, to use a (temporary) Boolean, in a numeric way.
(Numeric representation of Boolean, is even different in C (+1) / FB (-1).)

I'm probably going to use IIf(), which seems to fit the bill.

[deltarho[1859]]

Now, you are in denial.

I think that both of us should stop trying to score points - it is counterproductive.

Since, it is 'bad practice' IMO, to use a (temporary) Boolean, in a numeric way.
(Numeric representation of Boolean, is even different in C (+1) / FB (-1).)

I don't see a problem - True is -1 in FB and we are compiling FB. However, I will go along with IIf.

I'll look into your 'hack'

OK

[fxm]

This does not have distribution issues.

Actually, it does.

From FB's manual we have

Code: Select all

Function rnd_range (first As Double, last As Double) As Double
  Function = Rnd * (last - first) + first
End Function

This function is coherent with the definition of the interval of variation of Rnd():
in the range [0, 1[
(0 included and 1 excluded)

0 <= Rnd() < 1
as (last - first) >= 0,
0 <= Rnd() * (last - first) < (last - first)
first <= Rnd() * (last - first) + first < last
first <= rnd_range() < last
rnd_range() values are in the range [first, last[

The lower limit of the variation interval is always included and the upper limit is excluded.
(this makes it possible to apply the union operator on adjacent variation intervals and to obtain an equi-distributed result over the set)

[deltarho[1859]]

@fxm

Code: Select all

first <= rnd_range() < last

Agreed.

However, the values are not uniformly distributed in [0,1). There is a central bias giving a less likely return value at the range limits. My 'hack' seems to correct this.

[fxm]

Code: Select all

Function rnd_range (first As Double, last As Double) As Double
  Function = Rnd*(last - first -(first>=0)) - .5 + first
End Function 

For I As Integer = 1 To 50
  Print rnd_range (0, 1)
Next I

I do not understand.
We find values < 0 and > 1

[MrSwiss]

I do not understand. We find values < 0 and > 1

Simple:
You are using 'integer' data-type (instead of Double).

Just have a look at my figures ...

[deltarho[1859]]

Actually, it fails with Doubles as well.

My 'hack' is definitely required for LongInt - FB's rnd_range does not work with integer ranges.

I now wish that I had introduced only LongInt in the opening post and kept away from Doubles.

They have to be treated very differently. With an integer range of, say, 1 to 10 inclusive we have 10 possible results. The formula for i to j is 'j - i + 1' results. However, with a floating point range of, say, 1 to 10 if we consider the intervals 1 to 2, 2 to 3,..., 9 to 10 we only have 9 intervals.

My 'hack' does not work with Doubles but there is still a distribution issue with FB's rnd_range.

[MrSwiss]

My 'hack' does not work with Doubles but there is still a distribution issue with FB's rnd_range.

OK, I think, that it's fixable.
Since, first and last (what???) seem to be "unwisely" chosen names,
I've used: lo_lim and hi_lim (to better reflect their purpose).
A range, to my mind, is a limited space ... including the limits.

The REFERENCE implementation should be:

Code: Select all

Function rnd_range(ByVal lo_lim As Double, ByVal hi_lim As Double) As Double
    ' assure: inputs are correctly ordered
    If hi_lim < lo_lim Then Swap hi_lim, lo_lim
    Return Int( Rnd() * ((hi_lim + 1) - lo_lim) ) + lo_lim
End Function

Obtained results are (-5.0 to +5.0):

Code: Select all

-5                       9092718
-4.9                     9092718
-4.800000000000001       9092718
-4.700000000000001       9092718
-4.600000000000001       9092718
-4.500000000000002       9092718
-4.400000000000002       9092913
-4.300000000000003       9092913
-4.200000000000003       9092913
-4.100000000000003       9092913
-4.000000000000004       9092913
-3.900000000000004       9092913
-3.800000000000003       9092913
-3.700000000000003       9092913
-3.600000000000003       9092913
-3.500000000000003       9092913
-3.400000000000003       9090598
-3.300000000000003       9090598
-3.200000000000003       9090598
-3.100000000000003       9090598
-3.000000000000003       9090598
-2.900000000000003       9090598
-2.800000000000003       9090598
-2.700000000000002       9090598
-2.600000000000002       9090598
-2.500000000000002       9090598
-2.400000000000002       9090675
-2.300000000000002       9090675
-2.200000000000002       9090675
-2.100000000000002       9090675
-2.000000000000002       9090675
-1.900000000000002       9090675
-1.800000000000002       9090675
-1.700000000000002       9090675
-1.600000000000001       9090675
-1.500000000000001       9090675
-1.400000000000001       9094381
-1.300000000000001       9094381
-1.200000000000001       9094381
-1.100000000000001       9094381
-1.000000000000001       9094381
-0.9000000000000009      9094381
-0.8000000000000009      9094381
-0.700000000000001       9094381
-0.600000000000001       9094381
-0.500000000000001       9094381
-0.400000000000001       9086637
-0.300000000000001       9086637
-0.200000000000001       9086637
-0.100000000000001       9086637
-1.02695629777827e-015   9086637
 0.09999999999999898     9086637
 0.199999999999999       9086637
 0.299999999999999       9086637
 0.399999999999999       9086637
 0.499999999999999       9086637
 0.599999999999999       9090449
 0.699999999999999       9090449
 0.7999999999999989      9090449
 0.8999999999999989      9090449
 0.9999999999999989      9090449
 1.099999999999999       9090449
 1.199999999999999       9090449
 1.299999999999999       9090449
 1.399999999999999       9090449
 1.499999999999999       9090449
 1.599999999999999       9090057
 1.7                     9090057
 1.8                     9090057
 1.9                     9090057
 2                       9090057
 2.1                     9090057
 2.2                     9090057
 2.3                     9090057
 2.4                     9090057
 2.5                     9090057
 2.6                     9092298
 2.7                     9092298
 2.8                     9092298
 2.9                     9092298
 3                       9092298
 3.100000000000001       9092298
 3.200000000000001       9092298
 3.300000000000001       9092298
 3.400000000000001       9092298
 3.500000000000001       9090802
 3.600000000000001       9090802
 3.700000000000001       9090802
 3.800000000000001       9090802
 3.900000000000001       9090802
 4.000000000000001       9090802
 4.100000000000001       9090802
 4.2                     9090802
 4.3                     9090802
 4.4                     9090802
 4.499999999999999       9090802
 4.599999999999999       9088472
 4.699999999999998       9088472
 4.799999999999998       9088472
 4.899999999999998       9088472
 4.999999999999997       9088472

to high: 0   to low: 0   count: 100000000
out of range: 0

press a key to EXIT ...

[fxm]

A range, to my mind, is a limited space ... including the limits.

But this does not make it possible to apply the union operator on adjacent variation intervals and to obtain an equi-distributed result over the set.

[MrSwiss]

WTF, is a Union Operator?
The distribution is over the wole range equal (now).

[deltarho[1859]]

@MrSwiss

You have never published your code for testing Doubles.

With lo_lim = -5.0 how can we get a return value of, say, 2.3 from using

Code: Select all

Return Int( Rnd() * ((hi_lim + 1) - lo_lim) ) + lo_lim

In fact, if lo_lim does not have a decimal part I cannot see how anything can be returned with a decimal part.

But this does not make it possible to apply the union operator on adjacent variation intervals and to obtain an equi-distributed result over the set.

I must confess that I do not understand a word of that.

[fxm]

For example, when the two limits are included:
Let 'E[-5, 0]' be a set of random variables in [-5, 0]
Let 'E[0, +5]' be a set of random variables in [0, +5]
Let 'U' be the UNION operator over the 2 sets of random variables.
'E[-5, 0] U E[0, +5]' is not representative of 'E[-5, +5]', because 'E[-5, 0] U E[0, +5]' has a distribution peak for the 0 value (a factor 2).

Whereas if for example the upper limit is excluded:
'E[-5, 0[ U E[0, +5[' is representative of 'E[-5, +5['.