# 2019q3 Homework4 (compute-pi)
contributed by < `93i7xo2` >
###### tags: `sysprog2019`
Source: [2019q3 G06: compute-pi](https://hackmd.io/5C935B-JTLGAvHn0Zv-r5A?view)
:::success
## 作業要求
- [x] 在 GitHub 上 fork [compute-pi](https://github.com/sysprog21/compute-pi)，嘗試重現實驗，包含分析輸出
* 注意: 要增加圓周率計算過程中迭代的數量，否則時間太短就看不出差異
- [x] 閱讀 [error rate 分析](https://hackmd.io/@Ryspon/Hyd8pQPqx)，並思考我們該留意 error rate
- [ ] 用 gnuplot 作為 `$ make gencsv` 以外的輸出機制，預期執行 `$ make plot` 後，可透過 gnuplot 產生前述多種實做的效能分析比較圖表
- [ ] 可善用 [rayatracing](/s/B1W5AWza) 提到的 OpenMP, software pipelining, 以及 loop unrolling 一類的技巧來加速程式運作
- [ ] 詳閱前方 "Wall-clock time vs. CPU time"，思考如何精準計算時間
- [x] 除了研究程式以外，請證明 [Leibniz formula for π](https://en.wikipedia.org/wiki/Leibniz_formula_for_%CF%80)，並且找出其他計算圓周率的方法
- [ ] 思考如何不依賴 float 和 double 一類的浮點數型態來求圓周率，如 [不用函數庫和亂數 寫程式求 pi 值的方法](https://www.ptt.cc/bbs/Gossiping/M.1555612071.A.14C.html)
- [ ] 留意精準度: [你不可不知道的 double 十件事](https://www.ptt.cc/bbs/b97902HW/M.1222953645.A.AE5.html)
- [ ] 參照 [Bailey–Borwein–Plouffe formula](https://en.wikipedia.org/wiki/Bailey%E2%80%93Borwein%E2%80%93Plouffe_formula)
- [ ] 將你的觀察、分析，以及各式效能改善過程，並善用 gnuplot 製圖，紀錄於「[2019q3 Homework4 (作業區)](https://hackmd.io/@sysprog/Bkb5cDMuH)」
---
延伸思考:
- clock(), clock_gettime() 的使用
- time(), gettimeofday() 的使用
- 為什麼 clock_gettime() 結果飄忽不定？
- 為什麼 time() 和 gettimeofday() 不適合拿來作效能評比之用？
:::
## 實驗環境
```bash!
~$ lscpu
Architecture: x86_64
CPU op-mode(s): 32-bit, 64-bit
Byte Order: Little Endian
Address sizes: 39 bits physical, 48 bits virtual
CPU(s): 8
On-line CPU(s) list: 0-7
Thread(s) per core: 2
Core(s) per socket: 4
Socket(s): 1
NUMA node(s): 1
Vendor ID: GenuineIntel
CPU family: 6
Model: 94
Model name: Intel(R) Core(TM) i7-6700 CPU @ 3.40GHz
Stepping: 3
CPU MHz: 802.054
CPU max MHz: 4000.0000
CPU min MHz: 800.0000
BogoMIPS: 6816.00
Virtualization: VT-x
L1d cache: 32K
L1i cache: 32K
L2 cache: 256K
L3 cache: 8192K
NUMA node0 CPU(s): 0-7
~$ ldd --version
ldd (Ubuntu GLIBC 2.29-3ubuntu1) 2.29
~$ lsb_release -a
No LSB modules are available.
Distributor ID: Ubuntu
Description: Ubuntu 19.04
Release: 19.04
Codename: disco
```
## 重現實驗
[Error Rate]() 提及 loop=25 的實驗數值很不穩定，於是用 loop=250 重現實驗結果。
在 Makefile 裡， `plot` 這一段 rule 依賴 `gencsv` ，直接執行 `make plot` 產生圖片的同時一併產生測試數據 `result_clock_gettime.csv` 。
`METHOD ?=` 則用來指定計算方法，值為 BASELINE/LEIBNIZ/EULER 。
### BASELINE
Proof: 黎曼和
```cpp=
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;
}
```
### LEIBNIZ
證明參考 [Leibniz's Formula for Pi](https://proofwiki.org/wiki/Leibniz%27s_Formula_for_Pi)。從下面的 _Madhava–Leibniz series_ 開始推導:
$$
arctan(1) = \dfrac{\pi}{4} = 1 - \dfrac{1}{3} + \dfrac{1}{5} - \dfrac{1}{7} +\ ...
$$
首先積分下列[數列](https://proofwiki.org/wiki/Leibniz%27s_Formula_for_Pi/Lemma)
$$
\dfrac{1}{1+t^2} = 1 - t^2 + t^4 - t^6 + t^8 + ... + t^{4n} - \frac{t^{4n+2}}{1+t^2}
$$
從 $0$ 積分到 $x$, $0\leq{x}\leq1$
$$
\int_{0}^{x}\dfrac{1}{1+t^2}dt=x-\frac{x^3}{5}+\frac{x^5}{5}-\frac{x^7}{7}+...+\frac{x^{4n+1}}{4n+1}-R_n(x)
\\where\ R_n(x)=\int_{0}^{x}\dfrac{t^{4n+2}}{1+t^2}dt
$$
先看 $R_n(x)$ ，因為 $1\leq1+t^2$，得到
$$
0\leq{R_n(x)}\leq\int_{0}^{x}t^{4n+2}dt=\frac{x^{4n+3}}{4n+3}
$$
又因為
$$
\frac{x^{4n+3}}{4n+3}\leq\frac{1}{4n+3}, 0\leq{x}\leq1
$$
依據夾擠定理（Squeeze Theorem)，當 $n\rightarrow\infty, \frac{1}{4n+3}\rightarrow0$ ，於是得出下列式子
$$
\int_{0}^{x}\dfrac{1}{1+t^2}dt = x-\frac{x^3}{3}+\frac{x^5}{5}-\frac{x^7}{7}+\ ...
$$
而且$\frac{d}{dx}arctan(x)=\frac{1}{1+t^2}$
$$
arctan(x) = x-\frac{x^3}{3}+\frac{x^5}{5}-\frac{x^7}{7}+\ ...
$$
$x$ 代入 1 即為欲求結果
```cpp=
double compute_pi_leibniz(size_t N)
{
double pi = 0.0;
for (size_t i = 0; i < N; i++) {
double tmp = (i & 1) ? (-1) : 1;
pi += tmp / (2 * i + 1);
}
return pi * 4.0;
}
```
### EULER
```cpp=
double compute_pi_euler(size_t N)
{
double pi = 0.0;
for (size_t i = 1; i < N; i++) {
pi += 1.0 / (i * i);
}
return sqrt(pi * 6);
}
```
## Time function
- clock()
- clock_gettime()
- time()
- gettimeofday()
## Reference
- [OpenMP](https://computing.llnl.gov/tutorials/openMP/)
- [Ting19970 共筆](https://hackmd.io/ajtKm-SKTsOmQsUNegfk8g?view)
- [uccuser159共筆](https://hackmd.io/@uccuser159/2019q3Homework4_computepi)
- [Latex語法](https://zh.wikibooks.org/zh-tw/LaTeX/%E6%95%B0%E5%AD%A6%E5%85%AC%E5%BC%8F)
- [簡單學 makefile：makefile 介紹與範例程式 ](https://mropengate.blogspot.com/2018/01/makefile.html)
Sign in with Wallet
Connect another wallet