2016q3 Homework1 (raytracing)

# 2016q3 Homework1 (raytracing) contributed by <`janetwei`> reviewed by <`ierosodin`> - 可以利用gnuplot畫出優化前優化後的差異 - 利用pthread實作raytracing的平行化 - SIMD來加速整個資料的處理 - 程式碼和結果不要用圖片表達 ## 開發環境先載入ㄧ些開發工具 ```shell $ sudo apt-get update $ sudo apt-get install graphviz $ sudo apt-get install imagemagick ``` ## Raytracing 光線追蹤 (Ray tracing) 是三維電腦圖形學中的描繪演算法 (rendering algorithm)，跟蹤從眼睛發出的光線，而非是光源發出的光線，因此與物理現象有些差距 ![](https://i.imgur.com/1s1LJv8.png) convert ppm to png ![](https://i.imgur.com/d7nLfmt.png) check圖片輸出是否符合預期以 `make PROFILE=1` 重新編譯程式碼，並且學習 gprof ![](https://i.imgur.com/KuxWWyH.png) 可以發現某些function被呼叫了千萬次可以試著從這方面去優化程式 ## graphviz 因為graphviz需要使用特定格式,所以需要安裝gprof2dot轉換格式但因為需要用到pip,所以先安裝pyphon ``` $sudo apt install python python-pip $sudo pip install --upgrade pip $sudo pip install gprof2dot ``` 將gprof所得出的結果用graphviz表示 ![](https://i.imgur.com/LXGmGTo.png) 這是用編譯器最佳化選項得出的結果 ![](https://i.imgur.com/kuABrjP.png) 時間直接變成0.759s 因為將-O0- 變成 -Ofast- ![](https://i.imgur.com/1q20hd6.png) 但是我們的目標則是手動優化程式的結果比電腦更好（目標好遠大阿QQ） ## dot_product優化(loop unrolling) 由gprof分析 dot_product被呼叫了將近七千萬次，因此第一個想從這裡下手 ``` double dot_product(const double *v1, const double *v2) { double dp = 0.0; for (int i = 0; i < 3; i++) dp += v1[i] * v2[i]; return dp; } ``` ![](https://i.imgur.com/PHq62DB.png) 將迴圈展開 ``` double dot_product(const double *v1, const double *v2) { double dp = 0.0; for (int i = 0; i < 3; i++) dp = v1[0]*v2[0]+v1[1]*v2[1]+v1[2]*v2[2]; return dp; } ``` 在這個做完之後，我才發現這個叫做loop unrolling ![](https://i.imgur.com/EZZDORM.png) 發現已經減少到6.62s了，所以將其他function也用loop unrolling ![](https://i.imgur.com/jH7h0t2.png) >> 避免用圖片來表達程式碼和輸出的文字訊息，日後才能方便搜尋和複製，敬請配合 [name=jserv] 結果已經減少到5.51s了 ![](https://i.imgur.com/vaA1pl7.png) ## static inline static只對定義(有寫出程式碼)有用 inline 對編譯器有提示的作用,將function called轉成以下的程式碼,利用inline將程式碼展開,少了程式呼叫的成本但是我們在程式一開始時設定為-O0,最佳化被關掉所以inline沒有作用,因此我們要強制inline作用(force inline) 將`$ static inline`變成`$ static__attribute__((always_inline))` ![](https://i.imgur.com/UPuPcqR.png) 已經將時間減少為2.638s ``` Each sample counts as 0.01 seconds. % cumulative self self total time seconds seconds calls s/call s/call name 30.18 0.54 0.54 13861875 0.00 0.00 raySphereIntersection 29.62 1.07 0.53 13861875 0.00 0.00 rayRectangularIntersection 11.74 1.28 0.21 2110576 0.00 0.00 compute_specular_diffuse 7.26 1.41 0.13 1048576 0.00 0.00 ray_color 6.15 1.52 0.11 2110576 0.00 0.00 localColor 4.47 1.60 0.08 1048576 0.00 0.00 rayConstruction 4.19 1.68 0.08 4620625 0.00 0.00 ray_hit_object 2.24 1.72 0.04 1241598 0.00 0.00 reflection 1.68 1.75 0.03 1 0.03 1.79 raytracing 1.12 1.77 0.02 1241598 0.00 0.00 protect_color_overflow 0.56 1.78 0.01 1241598 0.00 0.00 refraction 0.56 1.79 0.01 1204003 0.00 0.00 idx_stack_push 0.28 1.79 0.01 113297 0.00 0.00 fresnel 0.00 1.79 0.00 2558386 0.00 0.00 idx_stack_empty 0.00 1.79 0.00 2520791 0.00 0.00 idx_stack_top 0.00 1.79 0.00 1048576 0.00 0.00 idx_stack_init 0.00 1.79 0.00 37595 0.00 0.00 idx_stack_pop 0.00 1.79 0.00 3 0.00 0.00 append_rectangular 0.00 1.79 0.00 3 0.00 0.00 append_sphere 0.00 1.79 0.00 2 0.00 0.00 append_light 0.00 1.79 0.00 1 0.00 0.00 calculateBasisVectors 0.00 1.79 0.00 1 0.00 0.00 delete_light_list 0.00 1.79 0.00 1 0.00 0.00 delete_rectangular_list ``` 可以知道接下來從raySphereIntersection著手優化 ## OpenMD 此網站上有著簡單介紹平行化處理的方法 https://kheresy.wordpress.com/2006/06/09/ 對於一般單一執行緒（single thread）的程式，多核心的處理器並沒有辦法提升它的處理效能；不過對於多執行緒（multi thread）的程式，就可以透過不同的核心同時計算，來達到加速的目的了！如果要在for迴圈平行化處理,在前面加上一行`$ #pragma omp parallel for`即可如果想利用 OpenMP 把 mian() 裡面的迴圈平行化處理呢？只需要修改成下面的樣子： ```C #include <omp.h> #include <stdio.h> #include <stdlib.h> void Test( int n ) { for( int i = 0; i < 10000; ++ i ) { //do nothing, just waste time } printf( "%d, ", n ); } int main(int argc, char* argv[]) { #pragma omp parallel for for( int i = 0; i < 10; ++ i ) Test( i ); system( "pause" ); } ``` 只要兩個或以上個for迴圈裏面的變數不是互相獨立,而且不是在平行化處理裡面宣告的變數要設定成private:在每個thread裡面都有個變數的副本在function raytracing中在3個for迴圈前加上對3個變數private 以及發現下面也可以做Loop unrolling ```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); #pragma omp parallel for num_threads(512)\ private(stk), private(d) \ private(object_color) 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 + 0 / factor, j * factor + 0 % 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; idx_stack_init(&stk); rayConstruction(d, u, v, w, i * factor + 1 / factor, j * factor + 1 % 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; idx_stack_init(&stk); rayConstruction(d, u, v, w, i * factor + 2 / factor, j * factor + 2 % 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; idx_stack_init(&stk); rayConstruction(d, u, v, w, i * factor + 3 / factor, j * factor + 3 % 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; //} } } } ``` 結果變成1.327s了 ``` janet@janet-MacBookAir:~/raytracing$ ./raytracing # Rendering scene Done! Execution time of raytracing() : 1.327288 sec ``` 因為執行`make PROFILE=1`會比直接`make`多花時間因此就試試看直接`make` ``` janet@janet-MacBookAir:~/raytracing$ ./raytracing # Rendering scene Done! Execution time of raytracing() : 0.788209 sec ``` 之後得到的時間為0.788s ## SIMD 可以同一時間點處理好幾筆資料

Syntax	Example	Reference
# Header	Header	基本排版
- Unordered List	Unordered List
1. Ordered List	Ordered List
- [ ] Todo List	Todo List
> Blockquote	Blockquote
Bold font	Bold font
Italics font	Italics font
~~Strikethrough~~	~~Strikethrough~~
19^th^	19^th
H~2~O	H₂O
++Inserted text++	Inserted text
==Marked text==	Marked text
[link text](https:// "title")	Link
![image alt](https:// "title")	Image
`Code`	`Code`	在筆記中貼入程式碼
```javascript var i = 0; ```	`var i = 0;`
:smile:		Emoji list
{%youtube youtube_id %}	Externals
$L^aT_eX$	L^aT_eX
:::info This is a alert area. :::	This is a alert area.