2016q3 Homework 1 (compute-pi)

contributed by <jeff60907>

tags: `jeff60907`

Reviewed by <`andy19950`>

有詳細介紹 user/sys/real time 的差異，I/O-bound CPU-bond 也有用簡單的中文解釋讓我們容易了解，但希望能夠把這些知識加入到效能分析內，比如 clock_gettime() 傳不同的參數分別拿到的是什麼樣的時間，用哪種時間算會比較精準。
有充分了解各種方法的執行成果，也可以嘗試看看用這些方法來證明看看別的演算法，像是萊布尼茲或蒙地卡羅之類的。
可以在 Makefile 修改迭代的次數，讓 gnuplot 畫出來的結果看起來比較好看，如果能取 95% 信賴區間更好。
最後是下面 compute_ci 的程式碼，在 hackMD 上看格式好像是錯的，我有嘗試修改但好像沒用，可能是在複製貼上的時候有貼到特殊字元吧。

開發環境

作業系統 Lubuntu 16.04
CPU:Intel(R) Core(TM) i7-3770 CPU @ 3.40GHz
L1d cache:             32K
L1i cache:             32K
L2 cache:              256K
L3 cache:              8192K

複習微積分概論
閱讀圓周率π 學長詳盡分析共筆
閱讀Function 其他演算法近期全部看一輪，惡補自己的未修過演算法的不足
2016q3 Homework 1 (compute-pi)

進行實驗之前

時間處理與 time 函式使用

What do 'real', 'user' and 'sys' mean

linux time 的指令

Real is wall clock time time from start to finish of the call. This is all elapsed time including time slices used by other processes and time the process spends blocked (for example if it is waiting for I/O to complete).

從一開始到結束所花費的時間，包含了 processes跟blocked 及 I/O的等待

User is the amount of CPU time spent in user-mode code (outside the kernel) within the process. This is only actual CPU time used in executing the process. Other processes and time the process spends blocked do not count towards this figure.

在user mode code 執行的CPU時間，只有CPU在程序的時間，其他processes跟blocked不計算再內。

Sys is the amount of CPU time spent in the kernel within the process. This means executing CPU time spent in system calls within the kernel, as opposed to library code, which is still running in user-space. Like 'user', this is only CPU time used by the process. See below for a brief description of kernel mode (also known as 'supervisor' mode) and the system call mechanism.

process 在kernel進行的CPU時間，像是 systme call

real time = sys + user ?

1.下指令 sleep（） real time > sys + user ，因為sleep 指令在 user mode 及 kernel mode 都沒有佔用 CPU
2.使用多核心計算可能會造成 real time < sys + user ，uesr是各個核心計算的總合

時間單位
- 秒用s表現，毫秒用ms,微秒用μs表示，納秒用ns表示，皮秒用ps表示
1. Wall-clock time 顧名思義就是掛鐘時間，也就是現實世界中實際經過的時間，是由 kernel 裡的 xtime 來紀錄，系統每次啟動時會先從設備上的 RTC 上讀入 xtime。
2. CPU time 指的是程序在 CPU 上面運行消耗 (佔用) 的時間，如果多核心運算 threads，CPU time 是「每條 thread 的使用時間加總」，單純計算 CPU time 是不夠的，因為執行時間可能還包括 I/O time、 communication channel delay、synchronization overhead…等等。

CPU bound 與 I/O bound

所謂 CPU bound，指的是程序需要大量的 CPU computation (對 CPU time 需求量大)，相比之下僅需少量的 I/O operation，效能取決於 CPU 速度。

所謂 I/O bound，需要大量 I/O operation (I/O request 的頻率很高或一次 I/O request 成本很高)，但僅需少量的 CPU computation，效能取決於 I/O device 速度。

很多時候我們可以使用 time 判斷程式是 CPU bound 還是 I/O bound。

real < user: 表示程式是 CPU bound，使用平行處理有得到好處，效能有提昇。
real ≒ user: 表示程式是 CPU bound，即使使用平行處理，效能也沒有明顯提昇，常見的原因有可能是其實計算量沒有很大，生成、處理、結束多執行緒的 overhead 吃掉所有平行處理的好處，或是程式相依性太高，不適合拿來作平行處理。
real > user: 表示程式是 I/O bound，成本大部分為 I/O操作，使得平行處理無法帶來顯著的效能提昇 ( 或沒有提昇)。

嘗試實驗

探討Makefile 及各檔案

閱讀Louei Lu

check: default
        time ./time_test_baseline
        time ./time_test_openmp_2
        time ./time_test_openmp_4
        time ./time_test_avx
        time ./time_test_avxunroll

$ make check 會執行default在各自執行這些指令
如果想單獨看某一檔案執行時間
$ time ./time_test_baseline

gencsv: default
        for i in `seq 100 5000 25000`; do \
                printf "%d," $$i;\
                ./benchmark_clock_gettime $$i; \
        done > result_clock_gettime.csv

可以看成for 初始100 每次 +5000 到 25000
然後把資料放進 .csv excel可讀檔

增加 make plot 到Makefile，先執行gencsv 在運作gnuplot

plot: gencsv
        gnuplot runtime.gp

要先把runtime.gp加進資料夾內，可參考phonebook

測試效能

$ make check

time ./time_test_baseline
N = 400000000 , pi = 3.141593
2.93user 0.00system 0:02.93elapsed 99%CPU (0avgtext+0avgdata 1788maxresident)k
0inputs+0outputs (0major+86minor)pagefaults 0swaps
time ./time_test_openmp_2
N = 400000000 , pi = 3.141593
2.93user 0.00system 0:01.47elapsed 198%CPU (0avgtext+0avgdata 1724maxresident)k
0inputs+0outputs (0major+86minor)pagefaults 0swaps
time ./time_test_openmp_4
N = 400000000 , pi = 3.141593
3.71user 0.00system 0:01.10elapsed 335%CPU (0avgtext+0avgdata 1604maxresident)k
0inputs+0outputs (0major+86minor)pagefaults 0swaps
time ./time_test_avx
N = 400000000 , pi = 3.141593
1.34user 0.00system 0:01.34elapsed 99%CPU (0avgtext+0avgdata 1696maxresident)k
0inputs+0outputs (0major+84minor)pagefaults 0swaps
time ./time_test_avxunroll
N = 400000000 , pi = 3.141593
1.18user 0.00system 0:01.18elapsed 99%CPU (0avgtext+0avgdata 1796maxresident)k
0inputs+0outputs (0major+84minor)pagefaults 0swaps

各自探討效能

TempoJiJi hackmd

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;
}

$ time ./time_test_baseline

N = 400000000 , pi = 3.141593

real	0m3.011s
user	0m3.008s
sys	0m0.000s

real ≒ user 表示程式是 CPU bound，即使使用平行處理，效能也不會明顯提昇。

openmp2

$ time ./time_test_openmp(2/4)















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;
}

N = 400000000 , pi = 3.141593

real	0m1.688s
user	0m3.368s
sys	0m0.000s

openmp4

N = 400000000 , pi = 3.141593

real	0m1.562s
user	0m4.748s
sys	0m0.000s

real < user: 表示程式是 CPU bound，使用平行處理有得到好處，效能有提昇。

avx

AVX 是一組 SIMD(Single Instruction Multiple Data)的指令集，擴展了原有的 MMX 及 Intel Streaming SIMD Extensions(Intel SSE)，增加了下列的 features

原來 128-bit 的 SIMD 暫存器擴展至 256-bit，未來可能會支援 512-bit 甚至 1024-bit
Three-operand，原來的 two-operand 運算是 A=A+B，將結果寫到來源暫存器，而新的是 A=B+C，舊的仍然存在
有一些指令可以支援四個暫存器的運算元，使得 code 更小，並且減少一些不必要的指令，提升速度
有一個新的 extension coding scheme(VEX)

$ time ./time_test_avx

N = 400000000 , pi = 3.141593

real	0m1.359s
user	0m1.360s
sys	0m0.000s

real ≒ user 表示程式是 CPU bound，即使使用平行處理，效能也不會明顯提昇。

avxunroll

$ time ./time_test_avxunroll

N = 400000000 , pi = 3.141593

real	0m1.246s
user	0m1.244s
sys	0m0.000s

real ≒ user 表示程式是 CPU bound，即使使用平行處理，效能也不會明顯提昇。

設定gnuplot格式進行圖表

reset    
                                                                    
set key width 2
set key left
set style fill solid

set yrange[0:]
set xrange[0:61000]
set grid

set title 'perfomance comparison'
set ylabel 'time(sec)'
set xlabel 'N'

set datafile separator "," 
set term png enhanced font 'Verdana,10'
set output 'runtime.png'


plot "result_clock_gettime.csv" using 1:2 with lines title 'baseline', \
'' using 1:3 with lines title 'omp_2', \
'' using 1:4 with lines title 'omp_4', \
'' using 1:5 with lines title 'avx', \
'' using 1:6 with lines title 'avxunroll'

runtime .gp 格式 xyrange可自行設定
圖一直編譯失敗，才發現CSV檔是以","做格子區分
set datafile separator "," 參考Vivan Lin

clokc_getime loop = 25 是25次取樣後才計算時間

之後將取樣loop 改為 1

這張圖發現 omp4 一開始抖動幅度很誇張，但是後面幾乎趨於平穩不在抖動

把疊帶次數縮短來看 omp4在500

發現把疊代數變多才看的出抖動的情形
這幾張baseline抖動的效果最大，omp4執行效能卻是最高的但是抖動程度也是很高
為什麼會造成這樣的抖動呢？ omp4其他同學效果都不是最好的 why?

實作學長分析共筆的error rate

先探討benchmark_clock_gettime.c內使用的 clock_gettime函式的原型為：

int clock_gettime(clockid_t clk_id, struct timespec *tp);

應用在計時的時候，第一個參數clk_id可填入CLOCK_REALTIME，而struct timespec* tp則是本函式回傳的結果。struct timespec的宣告如下：

struct timespec {
        time_t   tv_sec;        /* seconds */
        long     tv_nsec;       /* nanoseconds */
};

使用這種時間函數的話，可能可以精確到奈秒等級。
如果想要進一步確定可容許的精確度，則可以使用clock_getres，其原型如下：

int clock_getres(clockid_t clk_id, struct timespec *res);

time() returns the wall-clock time from the OS, with precision in seconds.
clock() seems to return the sum of user and system time.
clock_gettime(CLOCK_MONOTONIC, …) provides nanosecond resolution.
gettimeofday() returns the wall-clock time with (nominally) µs precision

參考C語言：使用clock_gettime計算程式碼的時間需求

Measure time in Linux - time vs clock vs getrusage vs clock_gettime vs gettimeofday vs timespec_get?

參考學長共筆的程式碼




double compute_error(double Pi){
    Pi=fabs(Pi-M_PI)/M_PI;
    return Pi;
}

error rate 分析

看到這張圖的時候以為自己是不是做錯了，怎麼avxunroll呈現鋸齒狀
然後9/29早上老師PO了一篇文原來也有同學發現N在1000結尾的時候會發生問題

把間距拉小看，可以發現AVX再約疊代一千次後幅度空間會變小
avxunroll 可以清楚的看到每約1000會有一個折反點，是什麼原因造成呢？

參考王紹華同學他提到avxunroll是因為avxunroll的 loop code執行一次等於 baseline 執行16次，所以會造成誤差，AVX 的版本是因為每個 for loop 只等同於 baseline 版本執行了4次，剛好可以整除資料點，才沒有算錯
ex：N=20，上述迴圈只會執行一次，而導致資料少了20 - 16 = 4筆

實作看看N疊代差異

seq 100 4800 500000

當設定間距為4800 / 16 的整數倍可以讓 avxunroll 不會因為誤差抖動，產生平滑曲線

seq 112 4800 500000

當起使點設定 112 /16 的整數倍卻讓avx 誤差完全消失，疑惑為什麼當起始點為16倍數反而造成這個現象？

信賴區間(Confidence Interval)

Steps:

計算平均值
計算標準差(standard deviation)
平均值的正負標準差乘以 2 就是 95% 信賴區間

標準差算法

加上信賴區間實作測

取樣 loop = 25

發現抖動情形沒有改變很多，在想程式碼是不是打錯
參考王紹華同學程式碼做修改
把取樣數 25個放進來計算 95%





































double compute_ci(double data[]) {
    double mean = 0;
    double SD = 0; // Standard Deviation
    double SDerror = 0; 
    double low,top;
    double result = 0;
    int count = 0;
  
    //compute mean
    for(int i = 0; i <loop; i++){
	mean += data[i];
    }
    
    mean = mean / loop;
	
    //compute Standard Deviation
    for(int i = 0; i < loop; i++) {
		SD += pow((data[i] - mean),2);
    }
    //Calculate standard error
    SDerror = sqrt(SD / (double)loop);

    //standard range
    low = mean - 2 * SDerror;
    top = mean + 2 * SDerror;

    // in  95%
    for(int i = 0; i < loop; i++) {
		if(data[i] >= low && data[i] <= top) {
	    	count++;
                result+= data[i];
		}
    }
    result = result / count;
    return result;
}

驗證 Leibniz formula for π








double compute_pi_leibniz(size_t N)
{
    double pi = 0.0;
    for (size_t i = 0; i < N; i++) {
        pi += pow(-1.0,i) / (2*i + 1);
    }
    return pi * 4.0;
}

測試了加上 -lm 讓math.h運作也無法，原來也有同學遇到這個問題
參考Vivan Lin
把平方拆成直接用

看起來leibniz 抖動效果看起來不怎麼明顯

未完其他之後補上

2016q3 Homework 1 (compute-pi)

tags: jeff60907

Reviewed by <andy19950>

開發環境

進行實驗之前

時間處理與 time 函式使用

嘗試實驗

探討Makefile 及各檔案

測試效能

各自探討效能

baseline

openmp2

openmp4

avx

avxunroll

實作學長分析共筆的error rate

error rate 分析

實作看看N疊代差異

信賴區間(Confidence Interval)

Steps:

標準差算法

加上信賴區間實作測

驗證 Leibniz formula for π

tags: `jeff60907`

Reviewed by <`andy19950`>