'''Leibniz formula implementation:''' 01 #include <iostream> 02 #include <iomanip> 03 04 #include <chrono> 05 06 // Function duplicated from [https://scs.senecac.on.ca/~gpu610/pages/workshops/w2.html Workshop 2 - BLAS] 07 void reportTime(const char* msg, std::chrono::steady_clock::duration span) { 08 auto ms = std::chrono::duration_cast<std::chrono::milliseconds>(span); 09 std::cout << msg << " - took - " 10 << ms.count() << " millisecs" << std::endl; 11 } 12 13 // Arctangent of 1; execute for ~n iterations to refine the result 14 long double arctan1( unsigned long long int it ) { 15 long double r = 0.0; 16 17 // v0.25; due to performing the operations in a different order, there are rounding issues... 18 for ( long double d = 0.0; d < it; d++ ) { 19 long double de = d * 4.0 + 1.0; 20 r += (1.0 / de) - (1.0 / (de + 2.0)); 21 } 22 23 return r; 24 } 25 26 int main( int argc, char* argv[] ) { 27 unsigned long long int n = std::stoull(argv[1], 0); 28 29 std::chrono::steady_clock::time_point ts, te; 30 ts = std::chrono::steady_clock::now(); 31 long double pi = 4.0 * arctan1( n ); 32 te = std::chrono::steady_clock::now(); 33 reportTime("Arctangent(1) ", te - ts); 34 35 // Maximum length of a long double is 64 digits; minus "3." gives 62 digits. 36 std::cout.precision(62); 37 std::cout << "Pi: " << std::fixed << pi << std::endl; 40 } Stage 1 - Big-O 

'''Monte-Carlo algorithm implementation:''' 01 #include <iostream> 02 #include <random> 03 04 #include <chrono> 05 06 // Function duplicated Duplicated from [https://scs.senecac.on.ca/~gpu610/pages/workshops/w2.html Workshop 2 - BLAS] 07 void reportTime(const char* msg, std::chrono::steady_clock::duration span) { 08 auto ms = std::chrono::duration_cast<std::chrono::milliseconds>(span); 09 std::cout << msg << " - took - " << 10 ms.count() << " millisecs" << std::endl; 11 } 12 13 int main(int argc, char* argv[]) { 14 std::chrono::steady_clock::time_point ts, te; 15 16 ts = std::chrono::steady_clock::now(); 17 unsigned long long int n = std::stoull(argv[1], 0), 18 totalCircle = 0; 19 20 int stride = 1000, 21 circleSize = n / stride; 22 23 unsigned int* circle = new unsigned int[circleSize]; 24 25 for (int i = 0; i < circleSize; i++) 26 circle[i] = 0; 27 28 std::random_device rd; 29 std::mt19937 mt(rd()); 30 std::uniform_real_distribution<long double> dist(0.0, 1.0); 31 te = std::chrono::steady_clock::now(); 32 reportTime("Init. ", te - ts); 33 34 ts = std::chrono::steady_clock::now(); 35 for (unsigned long long int i = 0; i < circleSize; i++) { 36 for (int j = 0; j < stride; j++) { 37 long double x = dist(mt), 38 y = dist(mt); 39 // if position is inside the circle... 40 if (x * x + y * y < 1.0) { 41 circle[i]++; 42 } 43 } 44 } 45 46 for (int i = 0; i < circleSize; i++) 47 totalCircle += circle[i]; 48 49 long double pi = 4.0 * ((long double) totalCircle) / ((long double) n); 50 te = std::chrono::steady_clock::now(); 51 reportTime("Drop points ", te - ts); 52 53 std::cout.precision(62); 54 std::cout << "Pi: " << std::fixed << pi << std::endl; 55 56 delete [] circle; 57 }