Generate Random Numbers Following a Normal Distribution in C/C++

Generate random numbers following a normal distribution in C/C++

There are many methods to generate Gaussian-distributed numbers from a regular RNG.

The Box-Muller transform is commonly used. It correctly produces values with a normal distribution. The math is easy. You generate two (uniform) random numbers, and by applying an formula to them, you get two normally distributed random numbers. Return one, and save the other for the next request for a random number.

C++ - generate random numbers following normal distribution within range

You don't specify the standard deviation. Assuming a standard deviation of 2000 for the given interval you can try this:

#include <iostream>
#include <random>

class Generator {
    std::default_random_engine generator;
    std::normal_distribution<double> distribution;
    double min;
    double max;
public:
    Generator(double mean, double stddev, double min, double max):
        distribution(mean, stddev), min(min), max(max)
    {}

    double operator ()() {
        while (true) {
            double number = this->distribution(generator);
            if (number >= this->min && number <= this->max)
                return number;
        }
    }
};

int main() {
    Generator g(7000.0, 2000.0, 1000.0, 11000.0);
    for (int i = 0; i < 10; i++)
        std::cout << g() << std::endl;
}

Possible output:

4520.53
6185.06
10224
7799.54
9765.6
7104.64
5191.71
10741.3
3679.14
5623.84

If you want to specify only the min and max values, then we can assume that the mean value is (min + max) / 2. Also we can assume that min and max are 3 standard deviations away from the mean value. With these settings we are going to throw away only 0.3% of the generated values. So you can add the following constructor:

Generator(double min, double max):
        distribution((min + max) / 2, (max - min) / 6), min(min), max(max)
    {}

and initialize the generator as:

Generator g(1000.0, 11000.0);

Random number from normal distribution in C++

Use a seed to initialize your generator. Here I am using a time-based seed.

#include <iostream>
#include <random>
#include <chrono>

using namespace std;

int main()
{
    unsigned seed = chrono::system_clock::now().time_since_epoch().count();
    default_random_engine generator(seed);
    normal_distribution<double> distribution(0.0, 1.0);

    cout << distribution(generator);
    return 0;
}

Random number within a range based on a normal distribution

A standard normal distribution has mean 0 and standard deviation of 1; if you want to make a distribution with mean m and deviation s, simply multiply by s and then add m. Since the normal distribution is theoretically infinite, you can't have a hard cap on your range e.g. (100 to 150) without explicitly rejecting numbers that fall outside of it, but with an appropriate choice of deviation you can be assured that (e.g.) 99% of your numbers will be within the range.

About 99.7% of a population is within +/- 3 standard deviations, so if you pick yours to be about (25/3), it should work well.

So you want something like: (normal * 8.333) + 125

How do I generate random numbers with a mean in C++?

Since you are flexible on distribution, an easy solution that still gives reasonable results, without having to do rejection logic, is a triangular distribution. I.e. you set the lower end of the triangle at 1,000, the upper end of the triangle at 20,000, and the tip of the triangle such that you get your desired mean of 9,000.

The wikipedia link above indicates that the mean of a triangular distribution is:

(a + b + c) / 3

where a and b are your lower and upper limits respectively, and c is the tip of your triangle. For your inputs, simple algebra indicates that c = 6,000 will give your desired mean of 9,000.

There is a distribution in C++'s <random> header called std::piecewise_linear_distribution that is ideal for setting up a triangular distribution. This needs only two straight lines. One easy way to construct such a triangular distribution is:

std::piecewise_linear_distribution<> dist({1000., 6000., 20000.},
                                          [](double x)
                                          {
                                              return x == 6000 ? 1. : 0.;
                                          });

Now you simply have to plug a URNG into this distribution and crank out the results. For sanity's sake it is also helpful to collect some statistics that are important according to your problem statement, such as minimum, maximum, and mean.

Here is a complete program that does this:

#include <algorithm>
#include <iostream>
#include <numeric>
#include <random>
#include <vector>

int
main()
{
    std::mt19937_64 eng;
    std::piecewise_linear_distribution<> dist({1000., 6000., 20000.},
                                              [](double x)
                                              {
                                                  return x == 6000 ? 1. : 0.;
                                              });
    std::vector<double> results;
    for (int i = 0; i < 100; ++i)
        results.push_back(dist(eng));
    auto avg = std::accumulate(results.begin(), results.end(), 0.) / results.size();
    auto minmax = std::minmax_element(results.begin(), results.end());
    std::cout << "size = " << results.size() << '\n';
    std::cout << "min = " << *minmax.first << '\n';
    std::cout << "avg = " << avg << '\n';
    std::cout << "max = " << *minmax.second << '\n';
}

This should portably output:

size = 100
min = 2353.05
avg = 8972.1
max = 18162.5

If you crank up the number of sampled values high enough, you will see convergence on your parameters:

size = 10000000
min = 1003.08
avg = 8998.91
max = 19995.5

Seed as desired.

Need a normally-distributed random number generator

The C++ standard library has a normal distribution class - just what you've asked for. Use it - and clip the value so it's between your minimum and maximum:

#include <algorithm> // for std::clamp()
#include <iostream>
#include <random>

int main() {
    std::random_device randomness_device{};
    std::mt19937 pseudorandom_generator{randomness_device()};

    auto mean = 1.5;
    auto std_dev = 0.5;
    auto min_allowed = 1.0;
    auto max_allowed = 2.0;
    std::normal_distribution<> distribution{mean, std_dev};
    auto sample = distribution(pseudorandom_generator);
    auto clamped = 
        // C++17 and later
        std::clamp(sample, min_allowed, max_allowed);
        // C++14 or earlier:
        // std::max(min_allowed,std::min(sample, max_allowed));
        // 

    std::cout 
        << "A value from a normal distribution with mean " << mean
        << " and standard deviation " << std_dev << ": "   << sample
        << "; when clamped to [" << min_allowed << ", " 
        << max_allowed << "], we get: " << clamped << "\n";
}

In terms of the distribution - this changes the measure so that the entire range (-infinity,1) is concentrated at 1, and similarly (2,infinity) is concentrated at 2. As commenters, suggest, there are other ways to interpret your request for an "approximately normal" distribution, such as re-sampling until you hit a value in the required range; or applying a continuous transformation which maps (infinity, infinity) to (1,2), e.g. x -> arctan(x). But you didn't specify what exactly you're after.

C++ fast normal random number generator

Most importantly, do you really need 100,000,000 random numbers simultaneously? The writing to and subsequent reading from RAM of all these data unavoidably requires significant time. If you only need the random numbers one at a time, you should avoid that.

Assuming that you do need all of these numbers in RAM, then you should first
profile your code if you really want to know where the CPU time is spent/lost.

Second, you should avoid unnecessary re-allocation and initialisation of the data. This is most easily done by using std::vector::reserve(final_size) in conjunction with std::vector::push_back().

Third, you could use a faster RNG than std::mt19937. That RNG is recommended when the quality of the numbers is of importance. The online documentation says that the lagged Fibonacci generator (implemented in std:: subtract_with_carry_engine) is fast, but it may not have a long enough recurrence period -- you must check this. Alternatively, you may want to use std::min_stdrand (which uses the linear congruential generator)

std::vector<double> make_normal_random(std::size_t number,
                                       std::uint_fast32_t seed)
{
  std::normal_distribution<double> normalDistr(0,1);    
  std::min_stdrand generator(seed);
  std::vector<double> randNums;
  randNums.reserve(number);
  while(number--)
    randNums.push_back(normalDistr(generator));
  return randNums;
}

---

Related Topics