# 2022q1 Homework2 (quiz2)
contributed by < `Uduru0522` >
###### tags: `linux2022`
## 問題一:Average Of Two Integer
最原始的寫法,逕直將兩樹相加除以二。倘若 `a + b > UINT_MAX` ,會造成 overflow。
```c
#include <stdint.h>
uint32_t average(uint32_t a, uint32_t b)
{
return (a + b) / 2;
}
```
> 不考慮以下所有效能增進,再者平均必定介於兩者之間,
> 其實也可以用一個 `uint64_t` 暫存 `a + b` 也可以避免 ovetflow 就是了(笑)
假設已知兩數的大小,可以改寫成下方的樣式,由於平均 `high` 及 `low` 必定介於兩者之間,可以避免 overflow;且已知大小可以避免在相減時造成 underflow。
```c
#include <stdint.h>
uint32_t average(uint32_t low, uint32_t high)
{
return low + (high - low) / 2;
}
```
### EXP1
首先,從嘗試從可免除 Overflow 的實作改寫成無除法的形式:
```c
#include <stdint.h>
uint32_t average(uint32_t a, uint32_t b)
{
return (a >> 1) + (b >> 1) + (EXP1);
}
```
從數式 $\text{average}=\text{low}+\frac{\text{high} - \text{low}}{2}$,可以推導出 $\text{average}=(\text{low} - \frac{\text{low}}{2}) + \frac{\text{high}}{2}=\color{red}{\frac{\text{low}}{2} + \frac{\text{high}}{2}}$
此時 `low` 及 `high` 之間不再進行計算,可以忽視他們的大小關係,可以 `a` 和 `b` 代替。
而 `a` 及 `b` 各自除以二則可以以 right-shift 代替除法。
然而,僅以 `(a >> 1) + (b >> 1)` 計算時由於 LSB 會被無條件捨棄,有可能造成 0.5 的誤差,當 `a` 和 `b` 雙方皆為奇數時,會導致計算結果相對於正確答案少了一的結果;因此,`EXP1` 便需要用來應對這個狀況。
:::info
舉例來說,當 `a = 3, b = 9`,`(a >> 1) + (b >> 1) = 1 + 4 = 5`,然而兩數平均皆為 6。
:::
綜上所述,`EXP1` 應該有以下的計算結果:
+ 當 `a`、`b` 皆為奇數結果為 **1**。
+ 其他狀況結果為 **0**。
可以得到 `(a & 1) & (b & 1)`,化簡可得 `a & b & 1`。
因此 `define EXP1 a & b & 1`。
### EXP2、EXP3
接著,嘗試更激進的改寫:
```c
uint32_t average(uint32_t a, uint32_t b)
{
return (EXP2) + ((EXP3) >> 1);
}
```
模仿 [binary adder](https://www.electronics-tutorials.ws/combination/comb_7.html),我們可以用 XOR 及 OR 運算替代掉加法:
$\require{cancel}\\
\begin{align}\\
\text{average}&=\frac{a+b}{2}\\
&= \frac{(a\ \mathtt{\text{XOR}}\ b)+(a\ \mathtt{\text{AND}}\ b) \ll 1}{2}\\
&= \frac{a\ \mathtt{\text{XOR}}\ b}{2}+\frac{(a\ \mathtt{\text{AND}}\ b)\cancel{\ll 1}}{\cancel{2}}\\
&= \frac{a\ \mathtt{\text{XOR}}\ b}{2}+(a\ \mathtt{\text{AND}}\ b)
\end{align}$
+ $(a\ \mathtt{\text{AND}}\ b)\times 2$ 計算**進位**,倘若在第 31 個 bit 的計算結果是 1,乘以2 (向左位移 1) 時便會造成溢位,因此利用取平均時除以二的動作消除這個位移。
+ $a\ \mathtt{\text{XOR}}\ b$ 計算**無進位加法**。
再將兩者分別除以 2 後相加,便可以得到兩數平均。接著探討這次相加有沒有可能造成溢位,以下以 $S$ 代稱 $a\ \mathtt{\text{XOR}}\ b\over 2$,$C$ 代稱 $(a\ \mathtt{\text{AND}}\ b)$:
+ ==溢位的原因為 $S$,$C$ 及進位的第 31 bit,三個之中有兩個以上是 1。==
- 例如:`b'01... + b'11...`
+ ==$S$ 由於除以 2 等同右移一位的效果,第 31 bit 必定是 0。==
==可以確定若有溢位,必定為 C 的第 31 bit 及 兩者第 30 bit 發生進位。==
從上述基礎開始推導:
+ 若 $C$ 第 31 bit 為 1,則 `a`,`b` 第 31 bit 皆為 1。
+ 由於 `a`,`b` 第 31 bit 皆為 1,$S$ 的第 30 bit 為 0。為使第 30 bit 發生進位,$C$ 的第 30 bit 以及第 29 bit 的進位必須為 1。
+ 若 $C$ 第 30 bit 為 1,則 `a`,`b` 第 30 bit 皆為 1。
+ 由於 `a`,`b` 第 30 bit 皆為 1,$S$ 的第 29 bit 為 0。為使第 29 bit 發生進位,$C$ 的第 29 bit 以及第 28 bit 的進位必須為 1。
+ ... ... ... ...
+ 若 $C$ 第 2 bit 為 1,則 `a`,`b` 第 2 bit 皆為 1。
+ 由於 `a`,`b` 第 1 bit 皆為 1,$S$ 的第 1 bit 為 0。為使第 1 bit 發生進位,$C$ 的第 1 bit 以及第 0 bit 的進位必須為 1。
+ 若 $C$ 第 1 bit 為 1,則 `a`,`b` 第 1 bit 皆為 1。
+ 由於 `a`,`b` 第 1 bit 皆為 1,$S$ 的第 0 bit 為 0,不符合前述條件,**可確定 $C$ 的第 31 bit 為 1 的狀況下第 30 bit 不可能進位**。
得證**本次加法不可能造成溢位**。
因此,`#define EXP2 a & b`,`#define EXP3 a ^ b`。
### 比較三者編譯後的輸出
編譯以下程式碼對三個函式利用 Compiler Explorer 進行比較:
+ 編譯器使用 x86-64 gcc 11.2。
+ 優化選項使用 `-O3`
+ 為了使三個函式的參數傳入規則相同,第一個版本加入了將變數更改為 `a > b` 的流程,以求公平。
```c
#include <stdint.h>
uint32_t average_1(uint32_t a, uint32_t b)
{
if(a < b){
a ^= b;
b ^= a;
a ^= b;
}
return b + (a - b) / 2;
}
uint32_t average_2(uint32_t a, uint32_t b)
{
return (a >> 1) + (b >> 1) + (a & b & 1);
}
uint32_t average_3(uint32_t a, uint32_t b)
{
return (a & b) + ((a ^ b) >> 1);
}
```
編譯輸出:
```asm
average_1:
cmp edi, esi
jnb .L2
mov eax, edi
mov edi, esi
mov esi, eax
.L2:
sub edi, esi
shr edi
lea eax, [rdi+rsi]
ret
average_2:
mov eax, edi
shr eax
mov edx, esi
shr edx
add eax, edx
and edi, esi
and edi, 1
add eax, edi
ret
average_3:
mov eax, edi
xor eax, esi
shr eax
and edi, esi
add eax, edi
ret
```
#### XOR swap / Normal swap
實驗途中,發現不論在 `-O1`、`-O2`、`-O3` 參數下,將 `a`、`b` 交換的流程都會將原始程式碼中的 XOR Swap 更改使用一個暫存器及 `mov` 指令來進行, 而不會用三個 `xor` 指令。
假設 `a` `b` 兩者原先都是存在於記憶體中,`mov` 指令的方式會以以下方式執行,==總共四個指令,四個記憶體操作==:
+ 讀取 `a` 到 `$1` (1)
+ 讀取 `b` 到 `$2`
+ 儲存 `$2` 到 `a`
+ 儲存 `$1` 到 `b`
而假設我們用 `xor` 指令進行處理,則會以以下方式執行,==總共七個指令,四個記憶體操作==:
+ 讀取 `a` 到 `$1`
+ 讀取 `b` 到 `$2`
+ **XOR** `a`、 `b` 到 `a`
+ **XOR** `a`、 `b` 到 `b`
+ **XOR** `a`、 `b` 到 `a`
+ 儲存 `$1` 到 `a`
+ 儲存 `$2` 到 `b`
但是此處的資料並不是從記憶體中提取,而是直接對 register 做操作,
編譯成 `xor %edi, %esi; xor %esi, %edi; xor %edi, %esi;` 貌似並不會影響執行速度,為何還需要更改成 `mov` 指令的樣式呢?
由於實在是無法整理出一個理由,決定直接到 StackOverflow 提問:
[Why is the XOR swap optimized into a normal swap using the MOV instruction?](https://stackoverflow.com/questions/71382441/why-is-the-xor-swap-optimized-into-a-normal-swap-using-the-mov-instruction)
而歸納出幾個原因:
+ **對 Instruction-level parallelism(ILP) 不友善:**
- 倘若我們使用 `mov` 指令,在將 `a`(`%edi`) 之內容移動至 `temp`(`%eax`) 之後,兩個 `mov` 指令可以平行執行;而 XOR swap 的方式則需要等待依序執行。
+ **Intel x86 架構下有更好的選擇:**
- *x86 指令集中存在 `XCHG`,作用為交換兩個 register 之中的內容*,速度比 XOR swap 的三個指令還要快。即使沒有多餘的 register 可以使用也應該使用這個指令。
+ **MOV-elimination 的加速:**
- intel 在 Ivy Bridge 以後的架構內建,在執行 `mov` 指令時與其直接移動 register 內的內容,而是修正內部的 register file 直接使其指向正確的值。比起執行速度很快的 `xor`,`mov` 直接**不執行**,速度更快。
此外,和 Clang 的輸出結果比較:
```asm
average_1: # @average_1
cmp edi, esi
mov eax, esi
cmovb eax, edi
cmovb edi, esi
sub edi, eax
shr edi
add eax, edi
ret
```
可以發現 CLang 利用 [Tail Duplication](http://elegantc.blogspot.com/2016/05/tail-duplication-with-llvm.html),有效的讓不論 `a > b` 或是 `a < b` 都可以在相同的指令數及 latencty 下執行,而 GCC 在 `-O1`, `-O2`, `-O3` 等配置下並未這麼做;因此我們可以手動進行這項優化:
```c
uint32_t average_1_taildup(uint32_t a, uint32_t b)
{
if(a < b){
return a + (b - a) / 2;
}
else{
return b + (a - b) / 2;
}
}
```
可以產出
```asm
average_1_taildup:
cmp edi, esi
jnb .L5
sub esi, edi
shr esi
lea eax, [rsi+rdi]
ret
.L5:
sub edi, esi
shr edi
lea eax, [rdi+rsi]
ret
```
犧牲些許的程式碼大小達成不論狀態皆為 5 個 uop 的結果。
:::info
根據 GCC 說明[(連結)](https://gcc.gnu.org/onlinedocs/gcc/Optimize-Options.html#:~:text=O2%2C%20-O3%2C%20-Os.-,-ftracer,-Perform%20tail%20duplication),Tail duplication 優化需要啟用 `-ftracer`,而它並不被包含在 `-O1`、`-O2`、`-O3`、`-Og`、`-Os` 的任意一個優化選項裡。添加之後則可以產出類似 `average_1_taildup()` 的編譯結果。
:::
### 研讀 linux/include/linux/average.h
```c=
#define DECLARE_EWMA(name, _precision, _weight_rcp) \
struct ewma_##name { \
unsigned long internal; \
}; \
static inline void ewma_##name##_init(struct ewma_##name *e) \
{ \
BUILD_BUG_ON(!__builtin_constant_p(_precision)); \
BUILD_BUG_ON(!__builtin_constant_p(_weight_rcp)); \
/* \
* Even if you want to feed it just 0/1 you should have \
* some bits for the non-fractional part... \
*/ \
BUILD_BUG_ON((_precision) > 30); \
BUILD_BUG_ON_NOT_POWER_OF_2(_weight_rcp); \
e->internal = 0; \
} \
static inline unsigned long \
ewma_##name##_read(struct ewma_##name *e) \
{ \
BUILD_BUG_ON(!__builtin_constant_p(_precision)); \
BUILD_BUG_ON(!__builtin_constant_p(_weight_rcp)); \
BUILD_BUG_ON((_precision) > 30); \
BUILD_BUG_ON_NOT_POWER_OF_2(_weight_rcp); \
return e->internal >> (_precision); \
} \
static inline void ewma_##name##_add(struct ewma_##name *e, \
unsigned long val) \
{ \
unsigned long internal = READ_ONCE(e->internal); \
unsigned long weight_rcp = ilog2(_weight_rcp); \
unsigned long precision = _precision; \
\
BUILD_BUG_ON(!__builtin_constant_p(_precision)); \
BUILD_BUG_ON(!__builtin_constant_p(_weight_rcp)); \
BUILD_BUG_ON((_precision) > 30); \
BUILD_BUG_ON_NOT_POWER_OF_2(_weight_rcp); \
\
WRITE_ONCE(e->internal, internal ? \
(((internal << weight_rcp) - internal) + \
(val << precision)) >> weight_rcp : \
(val << precision)); \
}
```
#### Exponentially weighted moving average (EWMA)
## 問題二
## 問題三
## 問題四
## 問題五
## 問題六
## 問題七:FizzBuzz
### 程式說明
#### 在不使用求餘運算下檢測數字可否被 3 或 5 整除
觀察程式碼,
可以注意到**以下的不等式用來得知一 32-bit 整數是否能被 d($d\in\{3,5\}$) 整除**:
> $S_{\mathtt{i}, 64}=\mathtt{\text{0xiiiiiiiiiiiiiiii}}, \mathtt{i}\in{[\text{0}_{16},\text{F}_{16}]}$
$n\times (S_{F,64}\div 3+1)\leq S_{F,64}\div 3$
經過化簡,可得
$n\times S_{5,64}+n\leq S_{5,64}$
接著,考慮三個不同 n 的值:整除、除三餘一、除三於二:
| $n\equiv$ | 0 (mod 3) | 1 (mod 3) | 2 (mod 3) |
| :-------: | :-------: | :-------: | :-------: |
| $n\times S_{F,64}+n$ | $\xrightarrow{\text{overflow}} n$ | $\xrightarrow{\text{overflow}}S_{5, 64}+n$ | $\xrightarrow{\text{overflow}}S_{A,64}+n$ |
利用 overflow 來取代除法運算中重複減去的過程,節省時間,
其中只有 $n\equiv 0(\text{mod 3})$ 時, 等式才成立。
:::info
此等式可以應用於檢測所有 $S_{F,64}$ 之因數 $k$ 且 $(S_{F,64}\div k) > L, L=$ **檢測範圍 bit 長度** 的
整數的因數上,舉例來說,
若要檢測是否可被 15 整除,$M_{15}$ 即可設成 $S_{1,64}+1$。
:::
**<div style="color:red">優點:$M_k$ 可提前計算,大幅減少除法運算的次數(僅一次)</div>**
#### 字串輸出優化
程式中字串輸出以『同個字串,不同起點和終點』來輸出,如以下表格:
| **字串內容** | F | i | z | z | B | u | z | z | %u | \0 |
| :-: | :-: | :-: | :-: | :-: | :-: | :-: | :-: | :-: | :-: | :-: |
| **起始點 (included)** | $0(mod 3)$ $0(mod 15)$| | | | $0 (mod 5)$ | | | | Others | |
| **終點 (excluded)** | | | | | $0(mod 3)$ | | | | $0 (mod 5)$ $0(mod 15)$ | Others |
#### 目標
由 `div5` 以及 `div3` 兩個值計算出**起點**以及**字串長**,`div5` $,$`div3` $\in\{0,1\}$
將上方表格轉換成數值:
| `div3` | `div5` | **起點** | `length` |
| :------: | :------: | :------: | :------: |
| 0 | 0 | 8 | $\geq$ 1 |
| 0 | 1 | 4 | 4 |
| 1 | 0 | 0 | 4 |
| 1 | 1 | 0 | 8 |
:::info
**起點方程式**
透過觀察,可以發現滿足以下方條件的方程式可以為解:
0. $\text{start_point}(0,0)=8$
1. **只要 `!div3` 便為0**
2. `div5` 會造成除以 2 的效果
可得 `start_point(div3, div5) = !div3 & (8 >> div5)`。(我的算法)
與題目中解答的相異之處為條件 1 的處理,題目使用的方式為
『若 `div3`,向右移 4-bit』,
無論運算值為何 (8 或 4) 皆會被消除。
#### `length` 方程式
透過觀察,可以發現滿足以下方條件的方程式可以為解:
0. $\text{length}(0,0)\geq 1$
1. `div3` 和 `div5` 任一為 1,結果為 4
2. `div3` 和 `div5` 同時成立會造成額外乘以 2 的效果
可得 `length = 2 << div3 << div5`。
:::