# 2016q3 Homework1 (raytracing)
contributed by <`abba123`>
## 開發環境
OS : ubuntu 16.04.1 LTS
除了 os 之外我們先把之後會用到的一些工具安裝起來吧
```
$ sudo apt-get update
$ sudo apt-get install graphviz
$ sudo apt-get install imagemagick
```
## 案例分析 Raytracing
### 目的
* 改寫檔案 math-toolkit.h 在內的函式實做，充分紀錄效能差異在共筆
* 可善用 POSIX Thread, OpenMP, software pipelining, 以及 loop unrolling 一類的技巧來加速程式運作
我們就先來實際跑一次吧
```
$ make
$ ./raytracing
```
耐心等一下後結果就出來啦
```
# Rendering scene
Done!
Execution time of raytracing() : 3.096112 sec
```
總共花了3.1秒的時間
這時會輸出一個 "out.ppm" 的圖檔，我們可以用 eog 開啟
```
$ eog out.ppm
```
### gprof分析
要使用 "gprof" 來分析，必須要在編譯檔案時，加上參數 '-pg'
我們來看一下我們的 Makefile
```c=
ifeq ($(strip $(PROFILE)),1)
PROF_FLAGS = -pg
CFLAGS += $(PROF_FLAGS)
LDFLAGS += $(PROF_FLAGS)
endif
```
這裡寫到假如我們在 make 後面加入參數 PROFILE=1，在編譯時就會加上 -pg
```
$ make PROFILE=1
```
先跑一次程式給 gprof 分析
```
$ ./raytracing
```
這時會產生一個 gmon.out 的文件作為分析用
開始分析吧
```
$ gprof -b raytracing gmon.out | less
```
這裡加了 less 是因為資訊太多，可以讓我們把訊息分很多頁，方便觀察
```
Flat profile:
Each sample counts as 0.01 seconds.
% cumulative self self total
time seconds seconds calls s/call s/call name
26.92 0.67 0.67 69646433 0.00 0.00 dot_product
14.87 1.04 0.37 56956357 0.00 0.00 subtract_vector
11.85 1.34 0.30 13861875 0.00 0.00 rayRectangularIntersection
10.45 1.60 0.26 31410180 0.00 0.00 multiply_vector
8.04 1.80 0.20 10598450 0.00 0.00 normalize
7.03 1.97 0.18 13861875 0.00 0.00 raySphereIntersection
4.42 2.08 0.11 17836094 0.00 0.00 add_vector
4.42 2.19 0.11 4620625 0.00 0.00 ray_hit_object
2.01 2.24 0.05 17821809 0.00 0.00 cross_product
2.01 2.29 0.05 1048576 0.00 0.00 ray_color
2.01 2.34 0.05 1 0.05 2.49 raytracing
1.61 2.38 0.04 2520791 0.00 0.00 idx_stack_top
1.21 2.41 0.03 4221152 0.00 0.00 multiply_vectors
0.80 2.43 0.02 2110576 0.00 0.00 compute_specular_diffuse
0.80 2.45 0.02 1048576 0.00 0.00 rayConstruction
0.40 2.46 0.01 2558386 0.00 0.00 idx_stack_empty
0.40 2.47 0.01 2110576 0.00 0.00 localColor
0.40 2.48 0.01 1204003 0.00 0.00 idx_stack_push
```
我們可以發現，時間大部分花在 dot_product() subtract_vector() rayRectangularIntersection()上，我們就針對他們做優化吧
## loop unrolling
為了減少執行時間，我們把 for 迴圈展開，減少 branch 產生
我們拿 subtract_vector() 當例子
```c=
static inline
void subtract_vector(const double *a, const double *b, double *out)
{
for (int i = 0; i < 3; i++)
out[i] = a[i] - b[i];
}
```
這邊我們已知 for 要做3次，我們便可以把它展開
```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];
}
```
雖然只降低一點點的時間，但由於此 function 呼叫次數多，加總起來也是很可觀
我們把可以做 unrolling 的迴圈都做一遍，然後再重新編譯執行一次
```
$ make
$ ./raytracing
```
```
# Rendering scene
Done!
Execution time of raytracing() : 2.291020 sec
```
結果降低了大約1秒左右，因為一個小小的改動，就讓時間大大降低了
## 強制inline
http://openhome.cc/Gossip/CppGossip/inlineFunction.html
在呼叫函式時會需要分配記憶空間因而需要額外的資源負擔
可以「建議」編譯器將之設定為「行內函式」（Inline function）
如果建議被採納，則該函式會自動在呼叫點展現為程式碼
雖然我們的 function 前面都有加 inline
但由於我們把編譯器最佳化關掉，所以此建議不會被採納
透過加入 '__attribute__ ((always_inline))' 來開啟
https://gcc.gnu.org/onlinedocs/gcc/Inline.html
```c=
__attribute__((always_inline)) 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];
}
```
```
$ make clean
$ make
$ ./raytracing
```
```
# Rendering scene
Done!
Execution time of raytracing() : 1.933347 sec
```
結果顯示，由於我們強至開啟 inline 時間又縮短了約0.3秒左右
## Openmp
接下來我們關注 raytracing() 這個 function
想辦法透過平行化，讓時間繼續縮短
```c=
void raytracing(uint8_t *pixels, color background_color,rectangular_node rectangulars, sphere_node spheres,light_node lights, const viewpoint *view,int width, int height)
{
point3 u, v, w, d;
color object_color = { 0.0, 0.0, 0.0 };
/* calculate u, v, w */
calculateBasisVectors(u, v, w, view);
idx_stack stk;
int factor = sqrt(SAMPLES);
for (int j = 0; j < height; j++) {
for (int i = 0; i < width; i++) {
double r = 0, g = 0, b = 0;
/* MSAA */
for (int s = 0; s < SAMPLES; s++) {
idx_stack_init(&stk);
rayConstruction(d, u, v, w,i * factor + s / factor,j * factor + s % factor,view,width * factor, height * factor);
if (ray_color(view->vrp, 0.0, d, &stk, rectangulars, spheres,lights, object_color,MAX_REFLECTION_BOUNCES)) {
r += object_color[0];
g += object_color[1];
b += object_color[2];
} else {
r += background_color[0];
g += background_color[1];
b += background_color[2];
}
pixels[((i + (j * width)) * 3) + 0] = r * 255 / SAMPLES;
pixels[((i + (j * width)) * 3) + 1] = g * 255 / SAMPLES;
pixels[((i + (j * width)) * 3) + 2] = b * 255 / SAMPLES;
}
}
}
}
```
我們發現，這裡面計算的東西互不影響，可以獨立運算
這也代表我們可以對這段程式做平行化
我們選用上手比較容易的Openmp來完成平行化
雖然容易上手，但缺有個缺點，就是不能處理巢狀的for
所以我們的任務是要把這三個for合併
http://stackoverflow.com/questions/28482833/understanding-the-collapse-clause-in-openmp
這篇下面有提到怎麼打散迴圈
```c=
double r[height][width],g[height][width],b[height][width];
#pragma omp parallel for
for(int n=0; n<height*width; n++) {
int i=n/width;
int j=n%height;
r[i][j]=0;
g[i][j]=0;
b[i][j]=0;
}
```
我們多建一個二維陣列來除存裡面的資訊
```c=
#pragma omp parallel for private(stk), private(d),private(object_color)
for (int n = 0; n < height*width*SAMPLES; n++) {
/* MSAA */
int i=n/(width*SAMPLES);
int j=(n/SAMPLES)%width;
int s=n%SAMPLES;
idx_stack_init(&stk);
rayConstruction(d, u, v, w, i * factor + s / factor,j * factor + s % factor,view,width * factor, height * factor);
if (ray_color(view->vrp, 0.0, d, &stk, rectangulars, spheres,lights, object_color,MAX_REFLECTION_BOUNCES)) {
r[i][j] += object_color[0];
g[i][j] += object_color[1];
b[i][j] += object_color[2];
} else {
r[i][j] += background_color[0];
g[i][j] += background_color[1];
b[i][j] += background_color[2];
}
pixels[((i + (j * width)) * 3) + 0] = r[i][j] * 255 / SAMPLES;
pixels[((i + (j * width)) * 3) + 1] = g[i][j] * 255 / SAMPLES;
pixels[((i + (j * width)) * 3) + 2] = b[i][j] * 255 / SAMPLES;
}
```
要使用 openmp 平行 for 迴圈只須在前面加上 '#pragma omp parallel for'
後面 private() 是為了確保每個 thread 裡的參數不會被影響到
在編譯時加上 '-fopenmp' 才能順利平行化
```
$ make
$ ./raytracing
```
```
# Rendering scene
Done!
Execution time of raytracing() : 1.006078 sec
Verified OK
```
透過平行化處理，又再次提昇了執行速度，1又快了秒左右。
## FUTURE WORK
* 用 POSIX Thread 自行切割
* 實做 AVX
Sign in with Wallet
Connect another wallet