tyrese
tyrese

Reputation: 11

Estimating the value of Pi using the Monte Carlo Method

I am writing a C program that will be able to accept an input value that dictates the number of iterations that will be used to estimate Pi.

For example, the number of points to be created as the number of iterations increases and the value of Pi also.

Here is the code I have so far:

#include <stdio.h>
#include <stdlib.h>
main()
{
    const double pp = (double)RAND_MAX * RAND_MAX;
    int innerPoint = 0, i, count;
    
    printf("Enter the number of points:");
    scanf("%d", &innerPoint);
    for (i = 0; i < count; ++i){
        float x = rand();
        float y = rand();
        
        if (x * x + y * y <= 1){
            ++innerPoint;
        }
        
        int ratio = 4 *(innerPoint/ i);
        
        printf("Pi value is:", ratio);
    }
}

Help fix my code as I'm facing program errors.

Upvotes: 1

Views: 810

Answers (2)

chux
chux

Reputation: 153303

rand() returns an integer [0...RAND_MAX].

So something like:

   float x = rand()*scale;  // Scale is about 1.0/RAND_MAX

The quality of the Monte Carlo method is dependent on a good random number generator. rand() may not be that good, but let us assume it is a fair random number generator for this purpose.

The range of [0...RAND_MAX] is RAND_MAX+1 different values that should be distributed evenly from [0.0...1.0].

((float) rand())/RAND_MAX biases the end points 0.0 and 1.0 giving them twice the weight of others.

Consider instead [0.5, 1.5, 2.5, ... RAND_MAX + 0.5]/(RAND_MAX + 1).

RAND_MAX may exceed the precision of float so converting rand() or RAND_MAX, both int, to float can incurring rounding and further disturb the Monte Carlo method. Consider double.

#define RAND_MAX_P1 ((double)RAND_MAX + 1.0)

// float x = rand();
double x = ((double) rand() + 0.5)/RAND_MAX_P1;

x * x + y * y can also incur excessive rounding. C has hypot(x,y) for a better precision sqrt(x*x + y*y). Yet here, with small count, it likely makes no observable difference.

//  if (x * x + y * y <= 1)
if (hypot(x, y <= 1.0))

Upvotes: 1

JANO
JANO

Reputation: 3066

I am sure it is not the best solution, but it should do the job and is similar to your code. Use a sample size of at least 10000 to get a value near PI.

As mentioned in the commenter: You should look at the data types of the return values functions give you.

#include <stdio.h>
#include <stdlib.h>
#include <time.h>

int main()
{
    
    // Initialize random number generation
    srand(time(NULL));
    
    int samples = 10000;
    int points_inside =0;

    // Read integer - sample size (use at least 10000 to get near PI)
    printf("Enter the number of points:");
    scanf("%d", &samples);
    
    
    
    for (int i = 0; i < samples; ++i){
        
        // Get two numbers between 0 and 1
        float x = (float) rand() / (float)RAND_MAX;
        float y = (float) rand() / (float)RAND_MAX;
        
        
        // Check if point is inside
        if (x * x + y * y <= 1){
            points_inside++;
        }
                
        // Calculate current ratio
        float current_ratio = 4 * ((float) points_inside / (float) i);
        
        printf("Current value of pi value is: %f \n", current_ratio);
    }
}

Upvotes: 0

Related Questions