Try   HackMD

2016q3 Homework1(Compute Pi)

tags: tundergod hw1 2016q3

contributed by <tundergod>

作業內容

作業要求

  • 在 GitHub 上 fork compute-pi，嘗試重現實驗，包含分析輸出
    • 注意: 要增加圓周率計算過程中迭代的數量，否則時間太短就看不出差異
  • 詳閱廖健富提供的詳盡分析，研究 rate rate 的計算，以及為何思考我們該留意
  • 用 phonebook 提到的 gnuplot 作為 $ make gencsv 以外的輸出機制，預期執行 $ make plot 後，可透過 gnuplot 產生前述多種實做的效能分析比較圖表
  • 可善用 rayatracing 提到的 OpenMP, software pipelining, 以及 loop unrolling 一類的技巧來加速程式運作
  • 詳閱前方 “Wall-clock time vs. CPU time”，思考如何精準計算時間
  • 除了研究程式以外，請證明 Leibniz formula for π，並且找出其他計算圓周率的方法
  • 將你的觀察、分析，以及各式效能改善過程，並善用 gnuplot 製圖，紀錄於「作業區」

開發環境

  • Linux 版本: Ubuntu 16.04 LTS

  • 硬體資訊：

Architecture:          x86_64
CPU op-mode(s):        32-bit, 64-bit
Byte Order:            Little Endian
CPU(s):                4
On-line CPU(s) list:   0-3
Thread(s) per core:    2
Core(s) per socket:    2
Socket(s):             1
NUMA node(s):          1
Vendor ID:             GenuineIntel
CPU family:            6
Model:                 69
Model name:            Intel(R) Core(TM) i5-4200U CPU @ 1.60GHz
Stepping:              1
CPU MHz:               1661.660
CPU max MHz:           2600.0000
CPU min MHz:           800.0000
BogoMIPS:              4589.54
Virtualization:        VT-x
L1d cache:             32K
L1i cache:             32K
L2 cache:              256K
L3 cache:              3072K
NUMA node0 CPU(s):     0-3

程式碼閱讀

  • computepi.c computepi.h是主程式
  • time_test.c 定義了N的大小
  • benchmark_clock_gettime利用lock_gettime()計算程式執行時間
  • make file中使用了-Dmacro依次執行不同的程式

原始數據

  • make check:
time ./time_test_baseline
N = 1000000000 , pi = 3.141593
11.01user 0.00system 0:11.01elapsed 99%CPU (0avgtext+0avgdata 1812maxresident)k
0inputs+0outputs (0major+85minor)pagefaults 0swaps

time ./time_test_openmp_2
N = 1000000000 , pi = 3.141593
12.22user 0.02system 0:06.14elapsed 199%CPU (0avgtext+0avgdata 1724maxresident)k
0inputs+0outputs (0major+89minor)pagefaults 0swaps

time ./time_test_openmp_4
N = 1000000000 , pi = 3.141593
24.18user 0.00system 0:06.13elapsed 394%CPU (0avgtext+0avgdata 1808maxresident)k
0inputs+0outputs (0major+93minor)pagefaults 0swaps

time ./time_test_avx
N = 1000000000 , pi = 3.141593
3.16user 0.00system 0:03.16elapsed 99%CPU (0avgtext+0avgdata 1932maxresident)k
0inputs+0outputs (0major+87minor)pagefaults 0swaps

time ./time_test_avxunroll
N = 1000000000 , pi = 3.141593
3.06user 0.00system 0:03.06elapsed 99%CPU (0avgtext+0avgdata 1784maxresident)k
0inputs+0outputs (0major+86minor)pagefaults 0swaps

不是很看得懂上面所顯示數據的後半部分，如果是手動輸入的時間只會有real time，user time,system time，如下：

time ./time_test_baseline 
N = 1000000000 , pi = 3.141593
real	0m11.085s
user	0m11.076s
sys	0m0.008s
  • 從以上的數據可以看到real time ≒ user time的時間差不多,該程式屬於CPU bound.所以測試中使用執行緒openmp時執行緒的數量越多執行時間反而越長,佔用的CPU也越高(2 therad:12.22; 4 thread:24.18).

計算PI的公式

實際操作

  • baseline:
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;
}
  • leibniz:
int main()
{
        double n = 500000000;
        double pi = 0.0;
        for(size_t i = 0; i < n; i++)
        {
                int sign = i % 2 == 0 ? 1 : -1;
                pi += (sign / (2.0 * (double)i + 1.0));

        }
        printf("pi = %lf\n",pi*4);
        return 0;

}
  • monte carlo:
int main()
{
   double n=500000000;
   double x,y;
   int count=0; /* # of points in the 1st quadrant of unit circle */
   double z;
   double pi;
   count=0;
   for(size_t i=0; i<n; i++) {
      x = (double)rand()/RAND_MAX;
      y = (double)rand()/RAND_MAX;
      z = x*x+y*y;
      if (z<=1) count++;
      }
   printf("pi = %lf \n",(double)count/n*4);
        return 0;
}

result of monte carlo: pi = 3.141539
  • monte carlo方法在計算5億個點後任有些許誤差,懷疑是rand函數的分布率非常差導致,會想辦法解決
    因爲c中的rand()函數會有固定輸出所以即使隨機多少個點都無法達到精準的準確度.
  • 對leibniz做avx 和 avx+unrolling;對monte carlo做openmp,結果:
time ./time_test_baseline
N = 500000000 , pi = 3.141593
5.76user 0.00system 0:05.76elapsed 99%CPU (0avgtext+0avgdata 1772maxresident)k
0inputs+0outputs (0major+84minor)pagefaults 0swaps

time ./time_test_openmp_2
N = 500000000 , pi = 3.141593
6.16user 0.02system 0:03.13elapsed 197%CPU (0avgtext+0avgdata 1824maxresident)k
0inputs+0outputs (0major+87minor)pagefaults 0swaps

time ./time_test_openmp_4
N = 500000000 , pi = 3.141593
11.43user 0.02system 0:03.16elapsed 362%CPU (0avgtext+0avgdata 1828maxresident)k
0inputs+0outputs (0major+94minor)pagefaults 0swaps

time ./time_test_avx
N = 500000000 , pi = 3.141593
1.77user 0.00system 0:01.77elapsed 99%CPU (0avgtext+0avgdata 1940maxresident)k
0inputs+0outputs (0major+90minor)pagefaults 0swaps

time ./time_test_avxunroll
N = 500000000 , pi = 3.141593
1.69user 0.00system 0:01.69elapsed 99%CPU (0avgtext+0avgdata 1724maxresident)k
0inputs+0outputs (0major+86minor)pagefaults 0swaps

time ./time_test_leibniz
N = 500000000 , pi = 3.141593
2.88user 0.00system 0:02.88elapsed 99%CPU (0avgtext+0avgdata 1724maxresident)k
0inputs+0outputs (0major+86minor)pagefaults 0swaps

time ./time_test_leibniz_avx
N = 500000000 , pi = 3.141593
1.85user 0.00system 0:01.85elapsed 99%CPU (0avgtext+0avgdata 1764maxresident)k
0inputs+0outputs (0major+86minor)pagefaults 0swaps

time ./time_test_leibniz_avx_unroll
N = 500000000 , pi = 3.141593
1.54user 0.00system 0:01.54elapsed 99%CPU (0avgtext+0avgdata 1720maxresident)k
0inputs+0outputs (0major+87minor)pagefaults 0swaps

time ./time_test_montecarlo
N = 500000000 , pi = 3.141539
15.24user 0.04system 0:15.29elapsed 99%CPU (0avgtext+0avgdata 1784maxresident)k
0inputs+0outputs (0major+84minor)pagefaults 0swaps

time ./time_test_montecarlo_openmp_2
N = 500000000 , pi = 3.141519
142.99user 186.48system 2:47.89elapsed 196%CPU (0avgtext+0avgdata 1816maxresident)k
0inputs+0outputs (0major+87minor)pagefaults 0swaps

time ./time_test_montecarlo_openmp_4
N = 500000000 , pi = 3.141506
129.12user 254.31system 2:07.48elapsed 300%CPU (0avgtext+0avgdata 1748maxresident)k
0inputs+0outputs (0major+92minor)pagefaults 0swaps

leibniz算式用avx和avx+unrolling做有明顯的效能增強,但是montecarlo利用openmp做卻變得非常緩慢,但是system time卻是最多的.

time ./time_test_montecarlo_openmp_4
N = 500000000 , pi = 3.141526
real	2m3.363s
user	2m1.964s
sys	3m58.884s

程式花費的時間分布:system time > real time = user time

   #pragma omp parallel num_threads(threads)
   {
       #pragma omp for private(x,y) reduction(+:z)
       for(size_t i=0; i<N; i++) {
          x = (double)rand()/RAND_MAX;
          y = (double)rand()/RAND_MAX;
          if (x*x+y*y <= 1) z++;
       }
   }

數據圖表1(baseline):

數據圖表2(leibniz):

數據圖表3(monte carlo method):

以上數據的sequences爲'seq 100 5000 100100'

Error Rate

  • references:
  • 爲什麼要考慮error rate?
    • 任何實驗都必須建立在正確或最小誤差的基礎上
    • monte carlo函式(如無法完全隨機分布)就是一個反面例子.
  • error rate的計算:
    • 利用定義在<math.h>中的M_PI與computepi中的計算結果對比
    • 但是在網絡上搜索pi的數值卻發現c所定義的pi精準度只達到小數15位
    ​pi found online = 3.14159265358979323846
​pi define in c  = 3.14159265358979311599
    • 誤差數據以compute_pi_baseline(),compute_pi_leibniz()和compute_pi_montecarlo()三個函數分別減去<math.c>裏面的M_PI測試.
    • 誤差數據圖表1:
    • 誤差數據圖表2:
    • 誤差數據圖表3:
    • 誤差數據圖表4:去掉logscale會發現

可以看到monte carlo method一直讀在跳動(隨機分布問題),而baseline和leibniz非常有趣(奇怪)的會以完全相同的誤差進行.並且精準值接近小數6位止(已以區近水平線的方式延展),在去掉logscale的第4張圖片會更加清楚的看到大約在N=6000的地方兩個演算法都以非常快的速度趨向最大精準值,而在6000之後對於pi精準小數的計算的效能非常差.
Lim Wen Sheng

Wall-clock time vs. CPU time

Measure time in Linux

  • time()
    • returns the wall-clock time from the OS, with precision in seconds.
    • 由於是從wall-clock time抓取的資料(1970開始經)所以在測試程式執行時間上如果有涉及NTP(Network Time Protocol),時區或是夏季節約時間,那就不是一個穩定的時間因爲會因爲而變成非單調遞增函數.

老師的公筆上說wall-clock time是從xtime上更新的,xtime以微妙microsecond(µs)計算,但是在一些網站看它是以秒來計算

  • clock()
    *returns an approximation of processor time used by the program.
    • 返回值类型是clock_t,它除以CLOCKS_PER_SEC(1000000)来得出时间,一般用两次clock函数来计算进程自身运行的时间.
    • 在32bit中clock()每72小時就會重新輪回一次(回傳同樣的值),232 = 72x60x100x100x100.
    • clock()並沒有考慮CPU被子進程使用的情況
    • 不能區分執行的權限是在user mode還是kernel mode
  • clock_gettime()
    • retrieve the time of the specified clock clk_id.
    • 精準度到達納秒,單調遞增,從RTC取值
    • 結構:
    ​struct timespec {
​time_t   tv_sec;        /* seconds */
​long     tv_nsec;       /* nanoseconds */
​};
  • getrusage()
    • 獲取進程所佔用系統的相關信息,獲取時間進準度可達微秒
    • 結構:
    ​struct rusage {
​​​​struct timeval ru_utime; /* user CPU time used */
​​​​struct timeval ru_stime; /* system CPU time used */
​​​​long   ru_maxrss;        /* maximum resident set size */
​​​​long   ru_ixrss;         /* integral shared memory size */
​​​​long   ru_idrss;         /* integral unshared data size */
​​​​long   ru_isrss;         /* integral unshared stack size */
​​​​long   ru_minflt;        /* page reclaims (soft page faults) */
​​​​long   ru_majflt;        /* page faults (hard page faults) */
​​​​long   ru_nswap;         /* swaps */
​​​​long   ru_inblock;       /* block input operations */
​​​​long   ru_oublock;       /* block output operations */
​​​​long   ru_msgsnd;        /* IPC messages sent */
​​​​long   ru_msgrcv;        /* IPC messages received */
​​​​long   ru_nsignals;      /* signals received */
​​​​long   ru_nvcsw;         /* voluntary context switches */
​​​​long   ru_nivcsw;        /* involuntary context switches */
​};
  • gettimeofday()
    • 回傳當前時間(有timeval和timezone兩個結構)
    • 回傳的時間會受到系統不連續跳躍的影響,所以並不適合用在計算執行時間上
  • mach_absolute_time()
    • 使用在OS X上的一個時間函數,精準度達納秒
  • 結論:
    • clock_gettime()幾乎可以在任何情況下使用,唯一的缺點是CPU在進行上下文切換時會發生微妙級別的小誤差.
    • time()和gettimeofday()因爲(有可能)不是單調遞增函數所以並不適合作benchmark

Xtime,Jiffies,Tick,Hz

  • xTime:從cmos電路上取得的時間,從1970年1月1日開始算到現在,Linux 核心也會決定每秒要發生幾次中斷 (透過常數 HZ) ，每次 timer interrupt，就會去更新 xtime。
  • Hz:linux kernel每固定的周期就會發出timer interrupt(防止某進程陷入無限迴圈中或是發生狀況讓CPU拿不回控制權)
  • Tick:
  • Jiffies:是記錄著從電腦開機到現在總共的時鐘中斷次數(tick).假如我們系統的頻率是200Mhz,那麼一次中斷的間隔是1s / 200,000，000Hz=0.000 000 005秒,所以理論上我們系統的精確度是5納秒。

延伸思考:

  • clock_gettime() 結果飄忽不定？ 
    ​The CLOCK_PROCESS_CPUTIME_ID and CLOCK_THREAD_CPUTIME_ID clocks are realized on
​many platforms using timers from the CPUs (TSC on i386, AR.ITC on Itanium). 
​These registers may differ between CPUs and as a consequence these clocks may
​return bogus results if a process is migrated to another CPU.
    • clock_gettime()能夠作用於電腦上所有的CPU,所以進程在CPU之間切換時會影響clock_gettime的準確度
    • 網絡上看到有人在dual intel5150上測試context switch對於不同的working set所花費的時間:

    • 解決方法:可以使用信賴區間解決數據飄忽不定的問題
      • 信賴區間
        • 信賴區間是統計學中對一個“概率樣本/參數”所進行的區間估計。信賴區間給出的是被測量參數的測量值的可信度。即，測量參數有一定的概率會落在該區間上。
      • 在對稱的常態分佈圖上，一般認爲機率樣本有68%會落在平均值加減變準差，而95%信賴區間則爲平均值加減2倍的標準差
      • 標準差公式：
      • 標準差的迷思：信賴區間所說的%並不是機率，而是一個信心值。並不代表隨機抽樣樣本一定有95%落在該區域上。(也稱爲68-95-99.7rule 或 empirical rule(經驗法則))

莱布尼茨公式的證明

Last changed by 
 
Lim Wen Sheng
0
1660
Published on HackMD