2016q3 Homework1 (Compute Pi)

tags: sysprog

開發環境

Ubuntu 14.04 LTS

$ lscpu

• CPU: Intel® Core™ i7-2600 CPU @ 3.40GHz × 8
• Mem: 8 GiB
• Cache:
  L1d 快取： 32K
  L1i 快取： 32K
  L2 快取： 256K
  L3 快取： 8192K

案例分析 - Compute Pi

Baseline

1. Formula of Pi

double compute_pi_baseline(size_t N)
{
    double pi = 0.0;
    double dt = 1.0 / N;                // dt = (b-a)/N, b = 1, a = 0
    for (size_t i = 0; i < N; i++) {
        double x = (double) i / N;      // x = ti = a+(b-a)*i/N = i/N
        pi += dt / (1.0 + x * x);       // integrate 1/(1+x^2), i = 0....N
    }
    return pi * 4.0;
}

time ./time_test_baseline
N = 400000000 , pi = 3.141593
38.15user 0.00system 0:26.11elapsed 146%CPU (0avgtext+0avgdata 1664maxresident)k
0inputs+0outputs (0major+91minor)pagefaults 0swaps

OpenMP

double compute_pi_openmp(size_t N, int threads)
{
    double pi = 0.0;
    double dt = 1.0 / N;
    double x;
    #pragma omp parallel num_threads(threads)
    {
        #pragma omp for private(x) reduction(+:pi)
        for (size_t i = 0; i < N; i++) {
            x = (double) i / N;
            pi += dt / (1.0 + x * x);
        }
    }
    return pi * 4.0;
}

2 threads

time ./time_test_openmp_2
N = 400000000 , pi = 3.141593
38.30user 0.00system 0:26.26elapsed 145%CPU (0avgtext+0avgdata 1664maxresident)k
0inputs+0outputs (0major+92minor)pagefaults 0swaps

4 threads

time ./time_test_openmp_4
N = 400000000 , pi = 3.141593
38.08user 0.00system 0:25.97elapsed 146%CPU (0avgtext+0avgdata 1628maxresident)k
0inputs+0outputs (0major+89minor)pagefaults 0swaps

Monte Carlo

隨機取點(x,y), 判斷是否再圓內

參考資料

double monte_carlo_pi(size_t N){
	double pi=0.0;
	const unsigned max = 10000;
	unsigned long i = 0;
	double x=0.0, y=0.0;
	unsigned long in_times=0;

	srand(time(NULL)); 
	for(i=0; i<N; i++){
		x=rand()%max;
		y=rand()%max;

		if(x*x + y*y < (double)(max*max))
			in_times++;
	}
	return (4.0 * (double)(in_times)/(double)(N));
}

time ./time_test_monte
N = 400000000 , pi = 3.141593
38.38user 0.00system 0:26.75elapsed 143%CPU (0avgtext+0avgdata 1664maxresident)k
0inputs+0outputs (0major+90minor)pagefaults 0swaps

Leibniz

Leibniz formula for π

double compute_pi_leibniz(size_t N)
{
    double pi = 0.0;
    for(int i=0; i<N; i++)
    {
        int temp = (i%2) ? -1 : 1;
        pi += (double) temp / (2*i+1);        
    }
    return pi * 4.0;
}

執行成果