http://www.fluentcpp.com/2017/03/02/the-pi-day-challenge-for-expressive-code-in-c/ | Estimate Pi using the Monte Carlo method with multi-threading

main.cpp

#include <iostream>
#include <random>
#include <vector>
#include <math.h>
#include <thread>
#include <mutex>
#include <iomanip>

#define PI 3.14159265359

double samplesInsideCircle = 0;
std::mutex samplesInsideCircle_mutex;

void findSamplesInsideCircle(int radiusPower, int samplePower, double &samplesInsideCircle)
{
	double radius = pow(10, radiusPower);
	double maxSamples = pow(10, samplePower);
	std::random_device random_device;
	std::mt19937 engine{ random_device() };
	std::uniform_real_distribution<double> distribution(0, radius);

	for (int i = 0; i < maxSamples; ++i)
	{
		double xDifference = distribution(engine) - radius;
		double yDifference = distribution(engine) - radius;
		double distance = sqrt(pow(xDifference, 2) + pow(yDifference, 2));

		if (distance < radius)
		{
			std::lock_guard<std::mutex> guard(samplesInsideCircle_mutex);
			samplesInsideCircle++;
		}
	}
}

int main(int argc, char** argv)
{
	std::vector<std::thread> threads;
	int radiusPower = 9;
	int samplesPower = 7;
	unsigned int maxThreads = std::thread::hardware_concurrency();

	for (unsigned int i = 0; i < maxThreads; i++)
	{
		threads.emplace_back(std::thread(findSamplesInsideCircle, radiusPower, samplesPower, std::ref(samplesInsideCircle)));
	}

	for (auto&& thread : threads)
	{
		thread.join();
	}
	
	double estimatedPi = (4 * samplesInsideCircle) / (pow(10, samplesPower) * maxThreads);
	std::cout << "Pi: " << estimatedPi << "\n";
	std::cout << "Difference between actual Pi: " << PI - estimatedPi << "\n";
	getchar();

	return 0;
}