ARM NEONにおける除算の精度と速度のトレードオフ / trade-off of  performance and accuracy of division in ARM NEON

#背景 / background
[前回の記事](https://qiita.com/sanmanyannyan/items/4d06b00078dd4abc4225)の続き。
今回は前回の補正方法を使った除算の速度計測。
In this article, I evaluated the performance of these division.

#性能測定 / performance evaluation
評価環境はraspberry pi 3 model b+を使った。
I used raspberry pi 3 model b+.
コンパイラはgcc 6.3.0を用いて計測した。
I used gcc 6.3.0.
詳細はコードみて。
If you want to go deep in details, Plz read the code.
[codes](https://github.com/sanmanyannyan/neon_div_performance_analysis)
一応比較のために加算/乗算も計測した。
I also measured the performance of add/mul as a reference.
それぞれの演算100,000回分の時間の計測を10,000回繰り返し、min, median, mean, maxをとった。
For these operation, I measured 100,000 times operation and iterate this for 10,000 times and get minimum, median, mean and max of execution time.

```
# execution results
iteration:100,000
trial:10,000
add
	 min:78.000000usec
	 med:79.000000usec
	mean:84.300500usec
	 max:356.000000usec
mul
	 min:78.000000usec
	 med:79.000000usec
	mean:78.766100usec
	 max:121.000000usec
div0: divisioin without correction
	 min:242.000000usec
	 med:243.000000usec
	mean:243.151000usec
	 max:285.000000usec
div2: division with 2nd order correction(1 step Newton-Raphson method)
	 min:821.000000usec
	 med:822.000000usec
	mean:822.104700usec
	 max:850.000000usec
div3: division with 3rd order correction
	 min:1121.000000usec
	 med:1122.000000usec
	mean:1122.309300usec
	 max:1150.000000usec
div4: division with 4th order correction(2 step Newton-Raphson method)
	 min:1350.000000usec
	 med:1350.000000usec
	mean:1351.053100usec
	 max:1556.000000usec
div8: division with 8th order correction(3 step Newton-Raphson method)
	 min:1864.000000usec
	 med:1865.000000usec
	mean:1865.662700usec
	 max:1909.000000usec
```

```div.s
# add part
.L30:
	add	r3, r3, #1
	cmp	r4, r3
	vadd.f32	q8, q8, q1
	vadd.f32	q7, q7, q3
	vadd.f32	q6, q6, q14
	vadd.f32	q5, q5, q12
	vadd.f32	q4, q4, q9
	vadd.f32	q8, q8, q2
	vadd.f32	q7, q7, q15
	vadd.f32	q6, q6, q13
	vadd.f32	q5, q5, q10
	vadd.f32	q4, q4, q0
	bne	.L30
```
```div.s
# mul part
.L34:
	add	r3, r3, #1
	cmp	r4, r3
	vmul.f32	q9, q9, q1
	vmul.f32	q7, q7, q3
	vmul.f32	q6, q6, q14
	vmul.f32	q5, q5, q11
	vmul.f32	q4, q4, q8
	vmul.f32	q9, q9, q2
	vmul.f32	q7, q7, q15
	vmul.f32	q6, q6, q12
	vmul.f32	q5, q5, q10
	vmul.f32	q4, q4, q0
	bne	.L34
```
```div.s
# div0 part
.L38:
	add	r3, r3, #1
	cmp	r4, r3
	vrecpe.f32	q10, q10
	vrecpe.f32	q9, q9
	vmul.f32	q10, q4, q10
	vmul.f32	q9, q4, q9
	vrecpe.f32	q10, q10
	vrecpe.f32	q9, q9
	vrecpe.f32	q5, q5
	vrecpe.f32	q7, q7
	vmul.f32	q5, q4, q5
	vmul.f32	q7, q4, q7
	vrecpe.f32	q5, q5
	vrecpe.f32	q7, q7
	vrecpe.f32	q6, q6
	vmul.f32	q10, q4, q10
	vmul.f32	q6, q4, q6
	vmul.f32	q9, q4, q9
	vrecpe.f32	q6, q6
	vmul.f32	q5, q4, q5
	vmul.f32	q7, q4, q7
	vmul.f32	q6, q4, q6
	bne	.L38
```

```div.s
# div2 part
.L42:
	add	r3, r3, #1
	cmp	r4, r3
	vrecpe.f32	q11, q9
	vrecpe.f32	q8, q7
	vrecps.f32	q9, q9, q11
	vrecps.f32	q7, q7, q8
	vmul.f32	q9, q9, q11
	vmul.f32	q8, q7, q8
	vmul.f32	q9, q10, q9
	vmul.f32	q8, q10, q8
	vrecpe.f32	q12, q9
	vrecpe.f32	q11, q8
	vrecps.f32	q9, q9, q12
	vrecps.f32	q8, q8, q11
	vmul.f32	q9, q9, q12
	vmul.f32	q11, q8, q11
	vmul.f32	q9, q10, q9
	vrecpe.f32	q8, q6
	vmul.f32	q7, q10, q11
	vrecps.f32	q6, q6, q8
	vrecpe.f32	q11, q5
	vmul.f32	q8, q6, q8
	vrecps.f32	q5, q5, q11
	vmul.f32	q8, q10, q8
	vmul.f32	q5, q5, q11
	vrecpe.f32	q11, q8
	vmul.f32	q5, q10, q5
	vrecps.f32	q8, q8, q11
	vrecpe.f32	q12, q5
	vmul.f32	q11, q8, q11
	vrecps.f32	q5, q5, q12
	vrecpe.f32	q8, q4
	vmul.f32	q5, q5, q12
	vrecps.f32	q4, q4, q8
	vmul.f32	q6, q10, q11
	vmul.f32	q8, q4, q8
	vmul.f32	q5, q10, q5
	vmul.f32	q8, q10, q8
	vrecpe.f32	q11, q8
	vrecps.f32	q8, q8, q11
	vmul.f32	q11, q8, q11
	vmul.f32	q4, q10, q11
	bne	.L42
```
#考察 / discussion
ラズパイ3のCPUがcortex-A53 @ 1.4GHzらしい。
cortex-A53の各命令のthroughput, latencyはわからなかったけど、[cortex-A55のよさげな資料](https://static.docs.arm.com/epm128372/20/arm_cortex_a55_software_optimization_guide_v2.pdf)は見つかったのでこれを参考にする。
この資料の32pから始まる表とその注釈を見てみると、ASIMDのFP加算命令は128bit幅でthroughput 1となっている。
加算命令100,000回分の時間は100,000[inst.] * 1[clock/inst.] / 1.4E+9[clock/sec] = 71usecで実測値とほぼ一致している。アセンブリを見ても加算命令以外ほとんど入っていないのでかなり正確に計測できてる。
前回の結果を鑑みると3次補正あたりがちょうどよさそう。
The machine has cortex-A53 @ 1.4GHz CPU.
I could not find throughput and latency of instructions in cortex-A53, but I could find [a useful document of cortex-A55](https://static.docs.arm.com/epm128372/20/arm_cortex_a55_software_optimization_guide_v2.pdf).
According to a table in 32p of the document, the throughput of FPP ADD instruction of ASIMD is 1. 
The execution time of 100,000 times Addition should be 100,000[inst.] * 1[clock/inst.] / 1.4E+9[clock/sec] = 71usec. This is close to the measured result. In the assembly code, there almost no not ADD instruction.  
We could measured accurate execution time.

#追記 / EDIT1
命令のthroughputとレイテンシーから考えるとdiv0, div2はadd/mulに比べて遅いのが気になったのでアセンブリを書き換えて性能を測ってみた。
div4, div8はめんどくさいから測ってない。

div0は乗算1回, 逆数推定命令1回でそれぞれのthroughputは1なので乗算の2倍程度なはず。ということでlatencyが隠れるように命令並び替えて計測したところ、無事に乗算の2倍ぐらいになって想定通り。

div2はニュートン法1ステップ分進める命令(VRECPS)を使っている。この命令のlatencyは8で5並列の演算じゃ隠し切れない。
div2は乗算2回, 逆数推定命令1回, VRECPS 1回なので乗算の4.5倍程度？
ということで測ってみたらまぁそのくらいだった。
VRECPSのlatencyがでかいのでかなりループアンロールしないと性能が出ないのか。

これくらいはgccが頑張って並び替えてくれると思ってたんだがそんな賢くないのね。

If we consider the throughput and latency, it is wired that div0 and div2 is quite slower than add/mul.
So I hand modified the assembly.

```
# execution results
iteration:100,000
trial:10,000
add
	 min:78.000000usec
	 med:79.000000usec
	mean:86.010500usec
	 max:338.000000usec
mul
	 min:78.000000usec
	 med:79.000000usec
	mean:78.761200usec
	 max:108.000000usec
div0
	 min:150.000000usec
	 med:150.000000usec
	mean:150.218200usec
	 max:177.000000usec
div2
	 min:350.000000usec
	 med:350.000000usec
	mean:350.349400usec
	 max:420.000000usec
```

```div_hand_opt.s
# div0 manually optimized
.L38:
	add	r3, r3, #1
	cmp	r4, r3
	vrecpe.f32	q5, q5
	vrecpe.f32	q6, q6
	vrecpe.f32	q7, q7
	vrecpe.f32	q9, q9
	vrecpe.f32	q10, q10
	vmul.f32	q5, q4, q5
	vmul.f32	q6, q4, q6
	vmul.f32	q7, q4, q7
	vmul.f32	q9, q4, q9
	vmul.f32	q10, q4, q10
	vrecpe.f32	q5, q5
	vrecpe.f32	q6, q6
	vrecpe.f32	q7, q7
	vrecpe.f32	q9, q9
	vrecpe.f32	q10, q10
	vmul.f32	q5, q4, q5
	vmul.f32	q6, q4, q6
	vmul.f32	q7, q4, q7
	vmul.f32	q9, q4, q9
	vmul.f32	q10, q4, q10
	bne	.L38
```
```div_hand_opt.s
# div2 manually optimized
.L42:
	add	r3, r3, #1
	cmp	r4, r3
	vrecpe.f32	q8, q4
	vrecpe.f32	q11, q5
	vrecpe.f32	q12, q6
	vrecpe.f32	q13, q7
	vrecpe.f32	q14, q9
	vrecps.f32	q4, q4, q8
	vrecps.f32	q5, q5, q11
	vrecps.f32	q6, q6, q12
	vrecps.f32	q7, q7, q13
	vrecps.f32	q9, q9, q14
	vmul.f32	q4, q4, q8
	vmul.f32	q5, q5, q11
	vmul.f32	q6, q6, q12
	vmul.f32	q7, q7, q13
	vmul.f32	q9, q9, q14
	vmul.f32	q4, q4, q10
	vmul.f32	q5, q5, q10
	vmul.f32	q6, q6, q10
	vmul.f32	q7, q7, q10
	vmul.f32	q9, q9, q10
	vrecpe.f32	q8, q4
	vrecpe.f32	q11, q5
	vrecpe.f32	q12, q6
	vrecpe.f32	q13, q7
	vrecpe.f32	q14, q9
	vrecps.f32	q4, q4, q8
	vrecps.f32	q5, q5, q11
	vrecps.f32	q6, q6, q12
	vrecps.f32	q7, q7, q13
	vrecps.f32	q9, q9, q14
	vmul.f32	q4, q4, q8
	vmul.f32	q5, q5, q11
	vmul.f32	q6, q6, q12
	vmul.f32	q7, q7, q13
	vmul.f32	q9, q9, q14
	vmul.f32	q4, q4, q10
	vmul.f32	q5, q5, q10
	vmul.f32	q6, q6, q10
	vmul.f32	q7, q7, q10
	vmul.f32	q9, q9, q10
	bne	.L42
```

#追記2
gccのオプションを試してみたところ-mtuneを指定しないと命令の並び替えをしてくれないらしい。
-mtune=cortex-a53を追加して試してみた。
I found mtune option of gcc.
I tried -mtune=cortex-a53.

```
iteration:100,000
trial:10,000
add
	 min:78.000000usec
	 med:79.000000usec
	mean:84.674700usec
	 max:288.000000usec
mul
	 min:78.000000usec
	 med:79.000000usec
	mean:79.066900usec
	 max:111.000000usec
div0
	 min:150.000000usec
	 med:151.000000usec
	mean:150.828100usec
	 max:237.000000usec
div2
	 min:350.000000usec
	 med:351.000000usec
	mean:351.729300usec
	 max:474.000000usec
div2 8
	 min:558.000000usec
	 med:559.000000usec
	mean:559.805000usec
	 max:687.000000usec
div3
	 min:515.000000usec
	 med:516.000000usec
	mean:516.764900usec
	 max:635.000000usec
div4
	 min:615.000000usec
	 med:617.000000usec
	mean:617.208100usec
	 max:746.000000usec
div8
	 min:852.000000usec
	 med:853.000000usec
	mean:853.973500usec
	 max:973.000000usec

```