RNG doubles, custom looking for feedback

[DarkLightHitomi]

So I am making a good generator class that uses a xorshift64* generation seeded by a xorshift+.

Well I made functions to return various number types and other suppirt functions.

One of the functions I just got to a working state is one that will return a double, including decimal values. So I'm curious what others think of it, and if there are any refinements I could make.

Note, No c&p available, so retyping this on my phone. Just in case you catch a mistyped something, let me know.

double rnD(double rng){//Function takes a double as a range maximum
bool one =FALSE;
If (rng ==1){one = TRUE; rng =2;}//A range of one breaks, so makes range two and sets flag for later correction
uint64_t x = rnRaw();//gets a raw value from the xorshift64* generator, rnRaw does seed checking so no need here
double y = 268435456;//for shifting the decimal point by 4 bytes
unint64_t d = rng*y;// shifts the rng into a whole number
x = x%d; //x becomes a value in range
double z = x/y; // shifts the decimal back by four bytes

if (one){ if (z>1){ z= z-1; //if rng had been 1 and corrected, undoes the correction 
}//close inner if scope 
}//close outer if scope 

return z;
}//close rnD() 

[Peter87]

Checking doubles for equality is often a mistake. What if someone pass in 0.99?

Is rnRaw() a random number? Why not simply divide the value from rnRaw() by the maximum value it can return and then multiply by rng to get the random number? I think that would give the randomness more precision.

[DarkLightHitomi]

Well, 0.99 doesnt cause a problem, but if you think about it, a number modulus 1 would always return the same result. Why that isn't solved by shifting the decimal is beyond me.

I also never thought of that divide by max method. I basically had to make it all from my own invention aside from the stripped down core of the xorshift generators.

Any idea what the max value of a uint64_t is?

[Peter87]

Well, 0.99 doesnt cause a problem, but if you think about it, a number modulus 1 would always return the same result. Why that isn't solved by shifting the decimal is beyond me.

Modulo 1 is usually not a problem, but modulo 0 is. If you pass in a value smaller than 1/y you have a problem.

Any idea what the max value of a uint64_t is?

std::numeric_limits<std::uint64_t>::max()

or

(std::uint64_t) -1

[DarkLightHitomi]

Hmm, tried your suggestion but now I get nothing but zeros.

using namespace std;

double rnD(double rng){
uint64_t x = rnRaw();
double z = x / (numeric_limits<uint64_t>::max());
return z*rng;

}//close rnD(double)

[DarkLightHitomi]

After some experimenting, It is the division that always results in a zero.

Printed max value so I also tried just dividing by the max value directly. Also tried dividing by max minus one, which also resulted in zeros. Tried removing just the last digit and got results but didn't remain in the range. Tried using max of double (since I was putting it into a double) but that gave massive wide values.

[DarkLightHitomi]

Now I explicitly casted the raw result from uint64_t into a double, then tried dividing by max size of double and multiplied by range, but get wildly out of range results which always seem to end with e+200 or more.

[DarkLightHitomi]

Peter87,
I think this trouble is an excellent reason to not do it your way.
:)

[Peter87]

If you divide two integers the result will be an integer. The decimal part is lost.

If you cast one of the operands to a double before doing the division it should work.

double z = (double) x / (numeric_limits<uint64_t>::max());

[DarkLightHitomi]

Okay, found a solution after some research,

double rnD(double rng){
uint64_t x = rnRaw();
double y = double(x) / double(numeric_limits<uint64_t>::max()/rng);
return y;
}//close function

Had to cast the terms of the expression. Also, moving range to the same expression to shorten the entire thing to three lines.

[DarkLightHitomi]

And of course, now you ninjad me. lol: )

Thanks for the info and help.

[JLBorges]

Consider making the custom generator a standard library compatible generator.
(minimum: meet uniform random number generator requirements,
if possible: meet random number engine requirements)

#include <cstdint>
#include <random>
#include <limits>

extern std::uint64_t rnRaw() ;

struct xorshift64_rng // meets uniform random number generator requirements
{
    using result_type = std::uint64_t ;
    static constexpr result_type min() { return 0 ; }
    static constexpr result_type max() { return std::numeric_limits<result_type>::max() ; }
    result_type operator()() const { return ::rnRaw() ; }
};

static xorshift64_rng my_rng ;

// generate a random double in the range [0,1)
double rnd()
{ return std::generate_canonical< double, std::numeric_limits<double>::digits >(my_rng) ; }
// http://en.cppreference.com/w/cpp/numeric/random/generate_canonical
// IS: note: Obtaining a value in this way can be a useful step in the process of transforming a value generated
//     by a uniform random number generator into a value that can be delivered by a random number distribution.
// IS: footnote:  b is introduced to avoid any attempt to produce more bits of randomness than can be held in RealType


// generate a random double in the range [lb,ub)
double rnd( double lb, double ub )
{ return std::uniform_real_distribution<double>(lb,ub)(my_rng) ; }
// http://en.cppreference.com/w/cpp/numeric/random/uniform_real_distribution

// generate a random double in the range [0,ub)
double rnd( double ub ) { return rnd( 0.0, ub ) ; }

[DarkLightHitomi]

Okay, you lost me. I'm not all that advanced in my knowledge (now feeling even more new than when taking the classes). Autism not helping either.

What do you mean by standard library compatible?

Why the struct?
Not really familiar with those so maybe I'm missing something about using a struct vs a class, but why do you have it all, what is the point behind "using return_type = uint64_t;" or the last line about "operator"?

What is with the extra stuff, like "generate_canonical" if the function can just return that anyway? (just have "return 1/ double(rnRaw());")

Isn't a high quality generator going to return a uniform real distribution? Wouldn't that be a major factor in determining quality?

[JLBorges]

> What do you mean by standard library compatible?

The standard library specifies the minimal requirements for a uniform random number generator. A generator which meets these requirements will interoperate in a well defined manner with the random number facilities provided by the library.

For instance, if my_random_number_generator is standard library compatible, the following would be possible:

#include <iostream>
#include <random>
#include <ctime>
#include <iomanip>
#include <algorithm>

// a standard library compatible uniform random number generator
struct my_random_number_generator // crude lcg: generates random numbers in the range [ 37, 547 ]
{
    my_random_number_generator( unsigned int seed ) : seed(seed) {}

    using result_type = unsigned int ;
    static constexpr result_type min() { return 37 ; }
    static constexpr result_type max() { return modulus + min() - 1 ; }
    result_type operator() () { return seed = ( multiplier * seed ) % modulus + min() ; }

    private:
        // note: only for exposition
        static constexpr unsigned int modulus = 511 ;
        static constexpr unsigned int multiplier = 48271 ;
        unsigned int seed = 0 ;
};

int main()
{
    my_random_number_generator rng( std::time(nullptr) ) ;
    std::cout << std::fixed << std::setprecision(1) ;

    // generate 10 integers uniformly distributed in [ -5000, +9000 ]
    std::uniform_int_distribution<int> distr( -5000, 9000 ) ;
    for( int i = 0 ; i < 10 ; ++i ) std::cout << distr(rng) << ' ' ;
    std::cout << '\n' ;

    // generate 10 floating point numbers uniformly distributed in [ 0, 1 )
    for( int i = 0 ; i < 10 ; ++i ) std::cout << std::generate_canonical<double,64>(rng) << ' ' ;
    std::cout << '\n' ;

    // generate 10 floating point numbers normal distribution with mean 13 and standard deviation 2.5
    std::normal_distribution<double> normal( 13.0, 2.5 ) ;
    for( int i = 0 ; i < 10 ; ++i ) std::cout << normal(rng) << ' ' ;
    std::cout << '\n' ;

    int array[] = { 0, 1, 2, 3, 4, 5, 6, 7, 8, 9 } ;
    // shuffle the contents of the array randomly
    std::shuffle( std::begin(array), std::end(array), rng ) ;
    for( int v : array ) std::cout << v << ' ' ;
    std::cout << '\n' ;
}

http://coliru.stacked-crooked.com/a/f78c5ad3123d94f7


> what is the point behind "using return_type = uint64_t;" or the last line about "operator"?

These are the requirements for a standard library compatible uniform random number generator: http://en.cppreference.com/w/cpp/concept/UniformRandomNumberGenerator


> Isn't a high quality generator going to return a uniform real distribution?

A uniform random number generator g of type G is a function object returning unsigned integer values such that each value in the range of possible results has (ideally) equal probability of being returned. [Note: The degree to which g’s results approximate the ideal is often determined statistically. —end note ] - IS

> What is with the extra stuff, like "generate_canonical" if the function can just return that anyway?
> (just have "return 1/ double(rnRaw());")

Let us say, the random number generator generates n-bit unsigned integers.
And the mantissa of double has m bits.

The rng can generate 2ⁿ distinct values.
To get a uniformly distributed double value in [0.0,1.0), we would need to pick one out of 2^m-1 possible distinct values. (m-1, assuming that we are not interested in denormal values.)

If n is less than m-1, to generate enough entropy, we would need to generate more than one n-bit random integer. std::generate_canonical does that in an optimal manner.

For instance, with our earlier toy my_random_number_generator:

#include <iostream>
#include <random>
#include <ctime>
#include <iomanip>
#include <limits>

// a standard library compatible uniform random number generator
struct my_random_number_generator // crude lcg: generates random numbers in the range [ 37, 547 ]
{
    my_random_number_generator( unsigned int seed ) : seed(seed) {}

    using result_type = unsigned int ;
    static constexpr result_type min() { return 37 ; }
    static constexpr result_type max() { return modulus + min() - 1 ; }
    result_type operator() () 
    { 
        static int n = 0 ;
        std::cout << ++n << ' ' ;
        return seed = ( multiplier * seed ) % modulus + min() ; 
    }

    private:
        // note: only for exposition
        static constexpr unsigned int modulus = 511 ;
        static constexpr unsigned int multiplier = 48271 ;
        unsigned int seed = 0 ;
};

int main()
{
    my_random_number_generator rng( std::time(nullptr) ) ;
    
    std::generate_canonical< float, std::numeric_limits<float>::digits >(rng) ; // calls rng() three times
    std::cout << '\n' ;

    std::generate_canonical< double, std::numeric_limits<double>::digits >(rng) ; // calls rng() seven times
    std::cout << '\n' ;

    std::generate_canonical< long double, std::numeric_limits<long double >::digits>(rng) ; // // calls rng() eight times
    std::cout << '\n' ;
}

http://coliru.stacked-crooked.com/a/cbd1553a5df4ee62

[DarkLightHitomi]

-Library compatibility, cool, but I'm making an entire class with all built in tools, so I'm not particularly worried about the library in this case.

-Distribution, It seems that many of the "better" generators return a uniform real distribution, including the ones I'm using.

-Canonical, according to my knowledge, 8.byte values are the largest integrals, so if someone is using a custom larger value, they'd have problems with this regardless and need to convert anyway (and can include methods for rng in that conversion), and so I'm not worried about it. Probably why the generators use uint64_t.

-For some odd reason (probably tiredness), I totally derped the op post with what generators I'm using. I'm building around a xorshift1024* seeded by a xorshift128+ not that it matters, but just if you're curious.

[JLBorges]

> It seems that many of the "better" generators return a uniform real distribution, including the ones I'm using.

3.2. GENERATING UNIFORM RANDOM NUMBERS

In this section we shall consider methods for generating a sequence of random fractions,
i.e., random real numbers U_n, uniformly distributed between zero and one. Since a computer can represent a real number with only finite accuracy, we shall actually be generating integers X_n, between zero and some number m; the fraction U_n = X_n/m will then lie between zero and one.

Usually m is the word size of the computer, so X_n may be regarded (conservatively) as the integer contents of a computer word with the radix point assumed at the extreme right, and U_n, may be regarded (liberally) as the contents of the same word with the radix point assumed at the extreme left.

Donald Knuth in 'Chapter 3 – Random Numbers, The Art of Computer Programming. Vol. 2: Seminumerical algorithms'

Donald Knuth 'Chapter 3 – Random Numbers, The Art of Computer Programming. Vol. 2: Seminumerical algorithms'
(Emphasis added)

[DarkLightHitomi]

Thanks for that. Don't know why in the blue blazes I was thinking 1/x would be uniform. Side effects of poor food and lack of sleep I guess.

Of coure, in the above example, m just needs to be the maximum possible result returned by the rng. I'm also guessing the radix point is the decimal point?

[JLBorges]

> I'm also guessing the radix point is the decimal point?

Binary point, in practice. http://en.cppreference.com/w/cpp/types/numeric_limits/radix

[DarkLightHitomi]

So decimal point, just slightly a misnomer as not all number systems are decimal, but really that point works exactly the same regardless of the number base.

Not really sure how the page you linked too was supposed to help though. Didn't seem to have anything to do with the previous text mentioning the radix point being on the left or right, while the page you just linked refers to number bases but no mention of a point between whole and partial numbers, much less placement of said point.