# 2017q Homework (raytracing)
contributed by <`tina0405`>,<`zhanyangch`>,<`LinRiver`>
### 開發環境
~~~
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: 58
Model name: Intel(R) Core(TM) i5-3230M CPU @ 2.60GHz
Stepping: 9
CPU MHz: 1200.024
CPU max MHz: 3200.0000
CPU min MHz: 1200.0000
BogoMIPS: 5188.41
Virtualization: VT-x
L1d cache: 32K
L1i cache: 32K
L2 cache: 256K
L3 cache: 3072K
NUMA node0 CPU(s): 0-3
~~~
### 以 `make PROFILE=1` 重新編譯程式碼，並且學習 [gprof](https://sourceware.org/binutils/docs/gprof/)
* 這邊使用 `make PROFILE=1` 是因為 Makefile 裡定義，而
`-pg` 的用法在下方GNU gprof 的原理提到
~~~=
ifeq ($(strip $(PROFILE)),1)
PROF_FLAGS = -pg
CFLAGS += $(PROF_FLAGS)
LDFLAGS += $(PROF_FLAGS)
endif
~~~
* GNU gprof 的原理，參考 [GNU gprof](https://sourceware.org/binutils/docs/gprof/Introduction.html#Introduction)
分析步驟 ：
1. 你應該編譯和連結你的程式使它可使用分析.往下解說如何編譯程式使之可分析
* `-pg` 的這個選項應該要跟著指令一起編譯和連結
* `-pg` 是編譯和連結的選項，如果沒使用就沒有 call-graph data
2. 我們應該要執行我們的程式使他產生分析的資料檔
* 一但程序被編譯成分析用，他還是會像原本那樣輸出，只是會多花時間蒐集和寫入資料
* 程序結束時會寫入 `gmon.out` 的檔案
3. 我們應該要使用 `gprof` 去分析資料
* 通常在編譯時加入 `-pg` ， 就代表著每一個函式會去呼叫 `mcount` (或 `_mcount`, 或 `__mcount`, 依靠作業系統和編譯器) 做為第一個操作之一
* 包含 profiling library 的 mcount 程序，通常利用檢查 stack frame 去找 child 的地址和返回 original parent 的地址
* mcount 是一個典型的簡短的組語 stub 程序，利用兩個參數 `frompc` 和 `selfpc` 呼叫 `__mcount_internal` ， `__mcount_internal` 是負責維護 in-memory call graph
記錄 frompc （從哪個 program counter 來） , selfpc （自己的 program counter）, 和調用次數。
* 首先必須先安裝 3 個套件 ， pip 是因為安裝 `gprof2dot` 會用到
~~~
$sudo apt install python python-pip
$sudo pip install --upgrade pip
$sudo pip install gprof2dot
~~~
#### 指令的使用
* 一開始如同上述說明，必須先執行一遍 `./raytracing` 將分析資訊儲存在 `gmon.out` 檔案中
* 有了分析資訊 `gmon.out` 就可以執行 gprof 分析 `gprof -p ./raytracing`
* -p 是只輸出函式調用圖
* -b 是不會輸出標題的描述
* -q 是只輸出時間消耗
~~~
Each sample counts as 0.01 seconds.
% cumulative self self total
time seconds seconds calls s/call s/call name
22.99 0.51 0.51 69646433 0.00 0.00 dot_product
17.58 0.90 0.39 56956357 0.00 0.00 subtract_vector
12.39 1.18 0.28 31410180 0.00 0.00 multiply_vector
9.47 1.39 0.21 10598450 0.00 0.00 normalize
9.01 1.59 0.20 13861875 0.00 0.00 rayRectangularIntersection
8.79 1.78 0.20 17836094 0.00 0.00 add_vector
4.96 1.89 0.11 4620625 0.00 0.00 ray_hit_object
4.51 1.99 0.10 13861875 0.00 0.00 raySphereIntersection
2.70 2.05 0.06 17821809 0.00 0.00 cross_product
1.80 2.09 0.04 2110576 0.00 0.00 compute_specular_diffuse
1.80 2.13 0.04 1048576 0.00 0.00 ray_color
1.13 2.16 0.03 3838091 0.00 0.00 length
0.90 2.18 0.02 2110576 0.00 0.00 localColor
0.68 2.19 0.02 4221152 0.00 0.00 multiply_vectors
0.45 2.20 0.01 2558386 0.00 0.00 idx_stack_empty
0.45 2.21 0.01 1241598 0.00 0.00 refraction
0.45 2.22 0.01 1048576 0.00 0.00 rayConstruction
0.00 2.22 0.00 2520791 0.00 0.00 idx_stack_top
0.00 2.22 0.00 1241598 0.00 0.00 protect_color_overflow
0.00 2.22 0.00 1241598 0.00 0.00 reflection
0.00 2.22 0.00 1204003 0.00 0.00 idx_stack_push
0.00 2.22 0.00 1048576 0.00 0.00 idx_stack_init
0.00 2.22 0.00 113297 0.00 0.00 fresnel
0.00 2.22 0.00 37595 0.00 0.00 idx_stack_pop
0.00 2.22 0.00 3 0.00 0.00 append_rectangular
0.00 2.22 0.00 3 0.00 0.00 append_sphere
0.00 2.22 0.00 2 0.00 0.00 append_light
0.00 2.22 0.00 1 0.00 0.00 calculateBasisVectors
0.00 2.22 0.00 1 0.00 0.00 delete_light_list
0.00 2.22 0.00 1 0.00 0.00 delete_rectangular_list
0.00 2.22 0.00 1 0.00 0.00 delete_sphere_list
0.00 2.22 0.00 1 0.00 0.00 diff_in_second
0.00 2.22 0.00 1 0.00 2.22 raytracing
0.00 2.22 0.00 1 0.00 0.00 write_to_ppm
~~~
* 文中提到 Graphviz 有 Dot 語言直連圖的能力，所以我們使用他將程式視覺化，以方便觀察需改善的函式
* 使用下列指令將程式視覺化
~~~
$ gprof ./raytracing | gprof2dot | dot -T png -o output.png
~~~
### 以 [gprof](https://sourceware.org/binutils/docs/gprof/) 指出效能瓶頸，並且著手改寫檔案 `math-toolkit.h` 的函式實做
* 執行時所花時間，有使用 `PROFILE=1`，跑 10 次的結果： 6.007 秒
~~~
Performance counter stats for 'sudo chrt -f 99 taskset -c 0 ./raytracing' (10 runs):
34 cache-misses # 27.528 % of all cache refs ( +- 7.22% )
125 cache-references ( +- 3.76% )
1676 instructions # 0.31 insn per cycle
5452 cycles ( +- 2.86% )
6.007348484 seconds time elapsed ( +- 0.14% )
~~~
* 觀看影片 [vedio1](https://www.youtube.com/watch?v=m1RmfOfSwno)，[vedio2](https://www.youtube.com/watch?v=V_rLyzecuaE)，參考 [Enhance raytracing program ： Team 4](https://embedded2016.hackpad.com/ep/pad/static/f5CCUGMQ4Kp) 的共筆
* 影片中提到如果我們開了最佳化（有使用 PROFILE=1）會得到以下結果，那表示我們目前未使用最佳化 6.007（sec），仍有許多方向要努力
~~~
tina@tina-X550VB:~/Hw/raytracing$ ./raytracing
# Rendering scene
Done!
Execution time of raytracing() : 0.741573 sec
~~~
* 列出程式熱區的前六名：由程式熱區當目標去改善，就算是少一行也有千萬級次數的影響
~~~
time seconds seconds calls s/call s/call name
22.99 0.51 0.51 69646433 0.00 0.00 dot_product
17.58 0.90 0.39 56956357 0.00 0.00 subtract_vector
12.39 1.18 0.28 31410180 0.00 0.00 multiply_vector
9.47 1.39 0.21 10598450 0.00 0.00 normalize
9.01 1.59 0.20 13861875 0.00 0.00 rayRectangularIntersection
8.79 1.78 0.20 17836094 0.00 0.00 add_vector
~~~
### Loop Unrolling [參考](http://book.51cto.com/art/200908/146356.htm)
Loop Unrolling 本身是一個相當基礎性的編譯器最佳化技術, 雖然他的方法看似直接了當, 但執行上具備相當大的成效. 其基本的觀念是把 for 迴圈中的表示式線性展開, 並且在每一次迭代中存取線性數據中的一小組來取代單獨一個. 如此一來程式在執行時可讓執行迭代次數減少.
以下為一針對彩色圖片取補色之程式來進行 Loop Unrolling, 其中 BYTE * pixel 為預處理之圖片像素數據, BYTE * tempPixel 為存放已處理之圖片, int width 和 height 為預處理圖片之寬與高
~~~
void Negative(BYTE * pixel, BYTE * tempPixel, int width, int height)
{
int i = 0;
int sum = width * height * 4;
memcpy(pixel, tempPixel, sum);
for(i = 0; i < sum; i++)
{
tempPixel[i] = 255 - tempPixel[i];
}
}
~~~
接著我們嘗試將迴圈次數變為原來三分之一, 因為每一像素具有 R, G, B 三分量, 透過以下改寫
~~~
void Negative(BYTE * pixel, BYTE * tempPixel, int width, int height)
{
int i = 0;
int sum = width * height * 4;
memcpy(pixel, tempPixel, sum);
for(i = 0; i < sum; i+=3)
{
tempPixel[i] = 255 - tempPixel[i];
tempPixel[i+1] = 255 - tempPixel[i+1];
tempPixel[i+2] = 255 - tempPixel[i+2];
}
}
~~~
然而需要特別注意的是, 根據預處理的數組的長度來選擇相對應適當的線性展開性. 例如增加暫存器儲存迭代的暫時變數, 因而導致效能降低. 除此之外, 更應考慮此展開性是否對時間與空間之開銷比例的影響
#### 使用共筆中 loop unrolling （有使用 PROFILE=1）
* 但 unrolling 也不是能隨意使用，因為一般現實狀況會有管線危害（pipeline hazard），要確定效能提升才可用 unrolling
* 函式次數最高者 `dot_product` 跑 10 次平均的結果 5.473391143 秒，比原版減少 0.5 秒
~~~c=
double dot_product(const double *v1, const double *v2)
{
double dp = (v1[0]*v2[0])+(v1[1]*v2[1])+(v1[2]*v2[2]);
return dp;
}
~~~
* 函式次數第二高者 `subtract_vector` 跑 10 次平均的結果 5.277608638 秒，再減少 0.2 秒
~~~c=
static inline
void subtract_vector(const double *a, const double *b, double *out)
{
out[0]=a[0]-b[0];
out[1]=a[1]-b[1];
out[2]=a[2]-b[2];
}
~~~
* Performance counter stats for 'sudo chrt -f 99 taskset -c 0 ./raytracing' (10 runs):
~~~
29 cache-misses # 24.138 % of all cache refs ( +- 8.68% )
119 cache-references ( +- 3.04% )
1676 instructions # 0.30 insn per cycle
5627 cycles ( +- 3.18% )
5.277608638 seconds time elapsed ( +- 0.44% )
~~~
* 函式次數第三高者 `multiply_vector` 跑 10 次的平均 5.277 秒，只減少 0.09
#### 使用共筆中 loop unrolling + macro 方法 （有使用 PROFILE=1）
* 改善呼叫函式次數最高者 `dot_product` ：
* `#define dot_product(a,b) ((a[0]*b[0])+(a[1]*b[1])+(a[2]*b[2]))`
* 執行跑 10 次的結果 ，執行時間： 4.176112849 秒 （比 `dot_product`+`subtract_vector` unrolling ）快了 1.1 秒
~~~
Performance counter stats for 'sudo chrt -f 99 taskset -c 0 ./raytracing' (10 runs):
51 cache-misses # 24.394 % of all cache refs ( +- 36.23% )
210 cache-references ( +- 44.10% )
3863 instructions # 0.38 insn per cycle ( +- 56.62% )
1,0112 cycles ( +- 45.37% )
4.176112849 seconds time elapsed ( +- 0.32% )
~~~
* 改善呼叫函式次數第二高者 `subtract_vector` ： 時間再加快 0.8 秒
~~~
Performance counter stats for 'sudo chrt -f 99 taskset -c 0 ./raytracing' (10 runs):
37 cache-misses # 29.299 % of all cache refs ( +- 10.34% )
126 cache-references ( +- 4.29% )
1676 instructions # 0.26 insn per cycle
6480 cycles ( +- 11.39% )
3.381460604 seconds time elapsed ( +- 0.68% )
~~~
* `dot_product` 和 `subtract_vector` 已經消失於圖中（函式）
* 改為 macro 的結果是時間減少但 cache-misses rate 增加，時間的提升原因來自於 69646433 次的`dot_product` 呼叫 + 56956357 次的 `subtract_vector` 呼叫
* 這個實驗中想要比較的是 macro 和 function 的差別
* function : 利用時間換取空間，大家共用一份程式區塊，因此所消耗的時間就是執行時期的 push 和 pop ，節省記憶體空間
* macro : 則是利用空間換取時間，前置處理器會將我們定義的 macro 替換各個區塊所以省下了 push 和 pop 的時間數，但造成空間的浪費
* 參考[前置處理器](http://epaper.gotop.com.tw/pdf/ael005600.pdf):
* C 語言的前置處理器是一個巨集處理器，主要用來處理C程式中含有 `#` 符號
1. 檔案含入功能： #include
2. 字串置換和巨集定義： #define、#undef
3. 條件編譯： #if...#elif...#else...#endif，#ifdef...#else...#endif，#ifndef...#else...#endif
* 動當編譯器在編譯程式之前前置處理器會自動啟：
* 可是即使時間從原本 6.007 變成 3.302 ， 加快 45% 但和最佳化後的 0.741 秒相差甚遠，因此再往 SIMD 和多執行緒的方面前進
### SIMD
* 參考 [春季班 Raytracing](https://hackmd.io/s/HkifAKiae#)
##### math-toolkit.h with SSE
* 解釋裡面常用指令，參考 [ Intel Intrinsics Guide ](https://software.intel.com/sites/landingpage/IntrinsicsGuide/#techs=SSE2&expand=3664) 這邊所使用的是 SSE2(一次可操作兩組雙精度符點數做運算)
* `__m128d _mm_loadu_pd (double const* mem_addr)`
* 操作方法 ：
~~~
dst [127:0] := MEM [mem_addr+127:mem_addr]
~~~
* 描述 : 載入 128 bits (由兩組雙精度符點數 64x2 = 128 bits) 從記憶體到 dst，mem_addr 不用對齊特定邊界
* `__m128d _mm_mul_pd (__m128d a, __m128d b)`
* 操作方法 ：
~~~
FOR j := 0 to 1
i := j*64
dst[i+63:i] := a[i+63:i] * b[i+63:i]
ENDFOR
~~~
* 描述 :如上述的 a 和 b 的 0~63 (64-bits) 相乘，放進 dst 0~63 (64-bits); a 和 b 的 64~127 (64-bits) 相乘，放進 dst 64~127 (64-bits)
* `__m128d _mm_add_pd (__m128d a, __m128d b)`
~~~
__m128d _mm_load_sd (double const* mem_addr)
__m128d _mm_add_pd (__m128d a, __m128d b)
__m128d _mm_div_pd (__m128d a, __m128d b)
void _mm_storeu_pd (double* mem_addr, __m128d a)
void _mm_store_pd1 (double* mem_addr, __m128d a)
__m128d _mm_sqrt_pd(__m128d a)
_mm_shuffle_pd (__m128d a, __m128d b, int imm8)
__m128d _mm_unpackhi_pd (__m128d a, __m128d b)
double _mm_cvtsd_f64(__m128d a)
~~~
* 先拿函式使用次數最高者 `dot_product` 來實測效能，所做的事其實就是 $(a[0] \times b[0])+(a[1] \times b[1])+(a[2] \times b[2])$
* force inline ： 在影片中提到我們關閉了最佳化 -O0 所以 inline 沒有效果，所以如果要強行執行的話，要加上 `__attribute__((always_inline))` 建議 CPU 執行它
* 所以用 SSE2 上述的指令集實作 `dot_product` 會像：
~~~
low high
取出 (v1[0] v1[1]) 打包放置 I0
取出 (v1[1] v1[2]) 打包放置 I1
取出 (v2[0] v2[1]) 打包放置 I2
取出 (v2[1] v2[2]) 打包放置 I3
|-----------|
I0*I2 (|v1[0]*v2[0]| v1[1]*v2[1]) 放置 T1
I1*I3 (|v1[1]*v2[1]| v1[2]*v2[2]) 放置 T2
T2high(|v1[2]*v2[2]| v1[2]*v2[2]) 放置 T3
|-----------|
框起來的部份全部低位相加是我們要的
缺點：會發現使用的空間和浪費的空間是 1 ： 1
~~~
* `dot_product` 程式碼部份：
~~~c=
__attribute__((always_inline)) static inline
double dot_product(const double *v1, const double *v2)
{
__m128d I0 = _mm_loadu_pd(v1);
__m128d I1 = _mm_loadu_pd(v1+1);
__m128d I2 = _mm_loadu_pd(v2);
__m128d I3 = _mm_loadu_pd(v2+1);
__m128d T1 = _mm_mul_pd(I0,I2);
__m128d T2 = _mm_mul_pd(I1,I3);
__m128d T3 = _mm_unpackhi_pd(T2,T2);
__m128d T4 = _mm_add_pd(T1,T2);
__m128d T5 = _mm_add_pd(T3,T4);
return _mm_cvtsd_f64(T5);
}
~~~
* 比較只改進 `dot_product` 部份，跑 10 次(都有 PROFILE=1)
* 原版 ： 6.263762944 seconds
* `dot_product`(unrolling + macro) : 4.176112849 seconds
* `dot_product`(sse) : 4.317184400 seconds
* 結果是比 unrolling+macro 慢，參考了春季班的共筆，提出可能性：
* 可能因為 128bit 運算 3 個 double ， 兩個兩個一組有所浪費
##### math-toolkit.h with AVX
* 參考 [ Intel Intrinsics Guide ](https://software.intel.com/sites/landingpage/IntrinsicsGuide/#techs=SSE2&expand=3664) 這邊所使用的是 AVX (一次可操作四組雙精度符點數做運算)
* 先拿函式使用次數最高者 `dot_product` 來實測效能：
~~~c=
static inline
double dot_product(const double *v1, const double *v2)
{
double temp[4] = {0, 0, 0, 0};
__m256d ymm0, ymm1;
ymm0 = _mm256_loadu_pd (v1); // ymm0 = {a0, a1, a2, 0}
ymm1 = _mm256_loadu_pd (v2); // ymm1 = {b0, b1, b2, 0}
ymm0 = _mm256_mul_pd (ymm0, ymm1); // ymm0 = {a0*b0, a1*b1, a2*b2, X}
const __m128d valupper = _mm256_extractf128_pd(ymm0, 1); // valupper = {a2*b2, 0}
ymm1 = _mm256_castpd128_pd256 (valupper); // ymm2 = {a2*b2, 0, X, X}
ymm0 = _mm256_hadd_pd (ymm0, ymm0); // ymm0 = {a0*b0+a1*b1, X, X, X}
ymm0 = _mm256_add_pd (ymm0, ymm1); // ymm0 = {a0*b0+a1*b1+a2*b2, 0, X, X}
_mm256_storeu_pd (temp, ymm0);
return *temp;
}
~~~
* 但其實這邊也可以看出來 AVX 沒寫好會比 SSE 更浪費空間
#### 整理 dot_product() 比較（有 make PROFILE=1，沒有 force inline）
* AVX 的效益比不上浪費空間所付出的代價
* 在 SIMD 和一般程式間轉換也需付出代價
### OpenMP
* 參考 [hugikun999 的共筆](https://hackmd.io/s/HyHhgcv6#)
以往在單核心處理器時期, 透過作業系統提供的 API 來建立執行緒. 然而在多核心的架構下, 為了讓各個處理器都能夠妥善被使用, 進而使用多執行緒程式設計來達到此目標. 透過 OpenMP 能夠取代使用作業系統之 API 來建立執行緒, 其中有以下三個特點:
* CPU 核心數量的擴展性
因為多核心程式設計時是需要考慮到程式效能會隨 CPU 數目增加而有所不同. 這其中要求程式中所建立的執行緒數量需要隨著 CPU 數目改變而變化. 然而使用 OpenMP 無須擔心此 CPU 數目擴展問題
* 便利性
使用 OpenMP 可將同一函式中的程式分成多個執行緒進行, 亦可將一 for 迴圈分解成多個執行緒進行
* 可移植性
各個主流的作業系統之執行緒 API 是互不相容, 然而 OpenMP 本身即是標準規範, 對於支持他的編譯器都是執行同一標準, 具備可移植性之特性
然而有幾項是使用 OpenMP 常見的錯誤使用:
* Use of locks without flush
* Use of ordered clause without ordered construct
* Try to change number of threads in parallel region after start of region
* Attempt to change loop variable while in #pragma omp for
* Use of critical when atomic would be sufficient
* Use of orphaned construct outside parallel region
* Use of unnecessary flush
* Use of unnecessary critical
### POSIX Threads
POSIX 執行緒, 亦可縮寫成 Pthreads, 為 [POSIX](https://zh.wikipedia.org/wiki/POSIX) 的執行緒標準, 定義了創建和操作執行緒的的一套 API, 一般用於 Unix-like POSIX 系統. 其 API 稱作 Pthreads, 大致共有100個函數調用且分為四類, 並以 'pthread_' 開頭:
* 執行緒管理
* 互斥鎖
* 條件變量
* 使用了互斥鎖的執行緒間的同步管理
### 目標
* [raytracing + SIMD]: 延續春季班 Raytracing 成果，確認 SIMD 和多執行緒得以運作
* [小組討論目標](https://hackmd.io/s/Hym6NzE-f#)
### 原作業要求
* 在 GitHub 上 fork [raytracing](https://github.com/sysprog21/raytracing)，並思考如何改善程式效能
* 以 `make PROFILE=1` 重新編譯程式碼，並且學習 [gprof](https://sourceware.org/binutils/docs/gprof/)
* 參考 [使用 GNU gprof 進行 Linux 平台的程式分析](http://os.51cto.com/art/200703/41426.htm)
* 以 [gprof](https://sourceware.org/binutils/docs/gprof/) 指出效能瓶頸，並且著手改寫檔案 `math-toolkit.h` 在內的函式實做，充分紀錄效能差異在共筆
* 注意: 請勿更動編譯器最佳化選項 `-O0` (關閉最佳化)
* 檔案 `math-toolkit.h` 定義的若干數學操作函式很重要，可參考 [Investigating SSE cross product performance](http://threadlocalmutex.com/?p=8)、[MathFu](https://google.github.io/mathfu/) 原始程式碼，以及 [2015q3 Homework #1 Ext](http://wiki.csie.ncku.edu.tw/embedded/2015q3h1ext)
* 可善用 POSIX Thread, OpenMP, software pipelining, 以及 loop unrolling 一類的技巧來加速程式運作
Sign in with Wallet
Connect another wallet