CODEWITHSUNDEEP
c++algorithmnumerical-methods
Ziezi
Ziezi

Reputation: 6467

Imprecision in Approximating Pi using Monte Carlo Method

Area of the circle = Pi * R^2 and the Area of the square = 4 * R^2. If we divide the area of the circle by the area of the square we get Pi / 4.

Let's have a square and an inscribed circle in it. We generate points with random coordinates and then count their number in each of the areas.
Then Pi = 4 * (# points in Circle) / (# points in Square).

Here is an attempt at approximating Pi with the above method: 

#include <iostream>     /* std::cout */
#include <iomanip>      /* std::setprecision */
#include <random>       /* std::uniform_int_distribution; std::mt19937 */

/* Check if point (x,y) is inside a circle with radius: r, at O(0,0). */
bool isInside (double x, double y, double r) { return (x*x + y*y) <= r*r; }

double approximatePi (double s, int sample_size)
{
    std::mt19937 gen;                   /* Generate random number in [-s/2 : s/2]. */
    std::uniform_int_distribution<double> dis(-s/2, s/2);

    int count = 0;                      /* Number of points in the circle. */ 
    for (int i = 0; i < sample_size; ++i)
    {
        double x = dis(gen); 
        double y = dis(gen);

        if (isInside(x, y, s/2))        /* Radius of inscribed circle = side / 2. */
        {
            ++count;
        }
    }
    /* Pi = 4 * (# points in Circle) / (# points in Square). */
    return (double) 4 * count / sample_size;    
}

int main()
{
    double side = 10.0;                 /* Square side. */
    int sample_size = 10000;            /* Number of tries. */

    std::cout <<"Pi ~ "<< std::fixed << std::setprecision(6) << approximatePi(side, sample_size) << '\n';
}

Expected result: Pi ~ 3.141592

Actual result: Pi ~ 2.611200

Why am I not getting the expected result? What am I doing wrong?

Upvotes: 1

Views: 97

Answers (1)

Drew Dormann
Drew Dormann

Reputation: 63946

The effect of std::uniform_int_distribution<double> is undefined behavior, because double is not an integer type.

Change it to std::uniform_real_distribution<double>.

Upvotes: 4

Related Questions