# 2016q3 Homework1 (compute-pi) #### contributed by <`andy19950`> ### 執行 make check ``` /*----- baseline -----*/ 16.23 user 0.00 system 0:16.23 elapsed 99% CPU /*----- OpenMP 2 threads -----*/ 16.21 user 0.00 system 0:08.11 elapsed 199% CPU /*----- OpenMP 4 threads -----*/ 17.70 user 0.00 system 0:04.75 elapsed 372% CPU /*----- AVX SIMD -----*/ 7.05 user 0.00 system 0:07.05 elapsed 99% CPU /*----- AVX SIMD + loop unrolling-----*/ 6.44 user 0.00 system 0:06.45 elapsed 99% CPU ``` - 可以發現使用 AVX SIMD 256bits暫存器來平行處理可以在執行時間上大幅縮短 --- ### 使用 gnuplot 繪製執行結果 - 首先需要新增.gp檔，我從第一個作業那邊複製過來做修改。 - 執行 make gencsv 會產生 result_clock_gettime.csv 裡面有執行各個funcion的時間 - 因為 gnuplot 不能用逗號 "，" 來做分隔所以必須修改 benchmark_clock_gettime.c - Makefile也要稍做修改，之後就可以使用 using 每筆資料都把第一列當作 x軸，輸出當 y軸。 - 如此執行就可以畫出每筆資料的 迭代次數(X) 以及 執行時間(Y) 的比較圖。 - logscale 可以對某個軸取 log ，我對 x 軸取 log10 讓圖看起來比較好看。 ``` clike= set logscale x 10 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' ``` --- ### 實做 Leibniz formula ![](https://i.imgur.com/neFZuxS.png) - 其實用上面的式子可以很簡單的實做出來 ```clike= double compute_pi_Leibniz(size_t N) { double pi = 0.0; double x; for (size_t n = 0; n < N; n++) { int temp = (n % 2) ? -1 : 1; x = (double) temp / (2 * n + 1); pi += x; } return pi * 4.0; } ``` --- ### 嘗試使用 AVX SIMD 來優化 Leibniz formula - AVX 256 bits 暫存器為 ymm0~ymm15。 - AVX 指令為 _mm256_op_suffix - op 我們有使用到 set, set1, add, mul, div, store - suffix 是參數的型態，因為我們用doulbe來計算所以使用 packed double(pd) - 因為使用 double 的關係 set 可以傳四個 doulbe (64bits) 參數把整個暫存器塞滿。 - 而 set1 則是傳一個 double 參數，暫存器的四個部份都儲存同一個參數。 ```clike= double compute_pi_Leibniz_avx(size_t N) { double pi = 0.0; register __m256d ymm0, ymm1, ymm2, ymm3, ymm4; ymm0 = _mm256_set_pd(1.0,-1.0,1.0,-1.0); ymm1 = _mm256_set1_pd(2.0); ymm2 = _mm256_set1_pd(1.0); ymm4 = _mm256_setzero_pd(); for (int i = 0; i <= N - 4; i += 4) { ymm3 = _mm256_set_pd(i, i+1.0, i+2.0, i+3.0); ymm3 = _mm256_mul_pd(ymm1,ymm3); ymm3 = _mm256_add_pd(ymm3,ymm2); ymm3 = _mm256_div_pd(ymm0,ymm3); ymm4 = _mm256_add_pd(ymm3,ymm4); } double tmp[4] __attribute__((aligned(32))); _mm256_store_pd(tmp, ymm4); pi += tmp[0] + tmp[1] + tmp[2] + tmp[3]; return pi * 4.0; } ``` - 因為要實現 (-1)^n^ 所以把 1 跟 -1 交錯存入ymm0 - 2*n + 1的部份用 ymm1, ymm2, ymm3來實做 - ymm4則當作答案的暫存器，最後使用 store 把 ymm4 存入記憶體 tmp 內 - 再把四個計算好的部份全部加起來 *4 就是答案 π !! --- ### 再使用 make check 跑跑看 ``` /*----- Leibniz -----*/ N = 1000000000 , pi = 3.1415926546 7.61 user 0.00 system 0:07.61 elapsed 99% CPU /*---- Leibniz + AVX SIMD -----*/ N = 1000000000 , pi = 3.1415926546 3.32 user 0.00 system 0:03.32 elapsed 100% CPU ``` - 結果證明萊布尼茲(Leibniz)的結論是對的。 - 而使用 AVX 來加速一樣有顯著的效果。 --- ### 使用 gnuplot 把所有結果畫成圖表 ![](https://i.imgur.com/8ywK7iQ.png) --- ### 其他估計 π 的方法 #### Nilakantha 級數 ----圖片來源：[wikiHow](http://zh.wikihow.com/%E8%AE%A1%E7%AE%97%E5%9C%86%E5%91%A8%E7%8E%87-Pi) ![](https://i.imgur.com/4BHA6Sk.png) - 據說使用電腦運算分析後， Nilakantha 級數會比 Leibniz formula 算的更快更準！ - 有時間我會寫程式來印證看看他的說法。 --- ### 結論 - 又把快遺忘的微積分拿出來複習一下，也發現很多方法可以來計算 π。 - gnuplot 又更會使用了，圖也畫的比較好看了。 - 使用電腦來做數學計算並驗證公式真的很有趣，也是我從來沒有想過的電腦使用方式。 - 本來前一個作業要用懂的 AVX SIMD 在這個作業終於搞懂他了。 --- ### Reference - [Intel Developer Zone](https://software.intel.com/en-us/node/524133) - [同學A共筆](https://embedded2015.hackpad.com/ep/pad/static/0aZx9N8YlMi) - [同學B共筆](https://hackmd.io/MYZgpgTAJmCsUFpYDNiwQFnMhBOAhmAGx5gYCMY+IUuRy5QA?view) - [WikiPedia](https://en.wikipedia.org/wiki/Leibniz_formula_for_%CF%80)