# 2024q1 Homework4 (quiz3+4)
contributed by < `david965154` >
### 第 3 周測驗題
#### 測驗 `1`
```c
int i_sqrt(int N)
{
int msb = (int) log2(N);
int a = 1 << (msb >> 1);
int result = 0;
while (a != 0) {
if ((result + a) * (result + a) <= N)
result += a;
a >>= 1;
}
return result;
}
```
第一題利用 2 的冪可以組成任何數的機制 (如同 4 bits 二進位 1111 可得出 0-31 之間的值,所做的僅需是計算每一位的值並做「總和」) ,而利用此方法可以大幅省去計算成本 $O(n)$->$O(logn)$ ,而實際上去執行時可以看出 msb 前面一半的值是完全不會用到的,也就是取 $N = 36$ 為例, a 會出現 32、16、8、4、2、1 ,不過 32、16、8 不可能會達成 `(result + a) * (result + a) <= N` 的條件,因此我將 msb 先右移 1 位再去做計算,並在每一次比較皆進行計算比較次數的行為,使比較次數減少約一半。
版本二
大致與版本一相同,僅差在不依賴函式 log2 並實作該功能。
版本三
直接跳到原式整理的地方
將 $N$ 使用二進位表示
$a_n$ 為二進位第 n 個 bit 代表的值
\begin{split} \\
P_m =& a_n+a_{n-1}+...+a_m \\
N^2 =& a_n^2+[2a_n+a_{n-1}]a_{n-1}+[2(a_n+a_{n-1})+a_{n-2}]a_{n-2}+...+[2(\displaystyle\sum_{i=1}^{n}a_i)+a_0]a_0 \rightarrow 1\\
=& a_n^2+[2P_n+a_{n-1}]a_{n-1}+[2P_{n-1}+a_{n-2}]a_{n-2}+...+[2P_1+a_0]a_0 \rightarrow 2\\
\end{split}
$P_0=a_n+a_{n-1}+...+a_0$ 即為 $N$
方法一、二即直接透過從 $a_n$ 試驗到 $a_0$ 湊出目標平方根 $N$ 二進位表示
試驗的過程為判斷 $P_m^2 \leq N^2$ 是否成立,而 $P_m = P_{m+1} + a_m$ ,因此會去看$P_m$ 加上 $a_m$ 是否達成條件,若成立 $P_m$ 即加上 $a_m$ ,反之則不加,以此逐漸逼近 $N$ 。
數學式表示如下:
\begin{cases}
P_m = P_{m+1} + 2^m, & \text{if $P_m^2 \leq N^2$}\\
P_m = P_{m+1}, & \text{otherwise}
\end{cases}
但是平方本身計算成本太高,因此往前一步所得去處理,不斷利用上個計算所得來取得這次計算所得,即遞迴式。
方法為此次差值 $X_m = N^2 - P_m^2$ 可改為上次差值 $X_{m+1}$ 減去 $P_m^2 - P_{m+1}^2$ (即 $P_m$ 之變化值,由是否加入 $a_m$ 而決定)
數學式表達如下:
- $X_m = N^2 - P_m^2 = X_{m+1} - Y_m$
- $Y_m = P_m^2 - P_{m+1}^2 = 2P_{m+1}a_m+a_m^2$
可以看到上數學式出現 $2P_{m+1}a_m+a_m^2$ 而 $a_m^2$ 與 2 的冪有關,可以作左右移達成計算,因此進一步作處理,將 $Y_m$ 拆解成 $c_m$ 與 $d_m$ 。
- $Y_m=\left.
\begin{cases}
c_m+d_m & \text{if } a_m=2^m \\
0 & \text{if } a_m=0
\end{cases}
\right.$
- $c_m = P_{m+1}2^{m+1}$
- $d_m = (2^m)^2$
再將 $c_m$ 與 $d_m$ 寫成遞迴的方式:
- $c_{m-1}=P_m2^m=(P_{m+1}+a_m)2^m=P_{m+1}2^m+a_m2^m=\begin{cases}
c_m/2+d_m & \text{if }a_m=2^m \\
c_m/2 & \text{if }a_m=0
\end{cases}$
- $d_{m-1}=\dfrac{d_m}{4}$
再來若要湊出 $N$ 的二進位表示法可用迴圈計算,而初始則先從:
1. $X_{n+1}=N^2$
:::warning
這裡應為 $N^2$ 而非原題目之 $N$
:::
第一步 $X_{n} = N^2 - P_n^2$ ,而 $P_n = a_n$ ,因此以 $X_{n+1}=N^2$ 表示。
2. $n$ 是最大的位數,使 $a_n^2=(2^n)^2= d_n^2 = 4^n \leq N^2$
即選取最合理之初始 n ,這裡將不會出現版本一、二會出現的 $a_n^2 > N^2$。
3. $c_n = 0$
從上面 $c_{m}$ 之遞迴式可以看到 $c_m=P_{m+1}2^{m+1}$ ,而在最初始時 $P_{n+1} = 0$ 因此 $c_n = 0$ 。
3. $d_n = 4^n$
從上面 $d_{m}$ 之遞迴式可以看到 $d_{m-1}=\dfrac{d_m}{4}$ ,經過 $n$ 次計算時 $d_{0} = 0$ ,因此原始 $d_{n} = 4^n$ 。
迴圈結束條件:
因為 $d_m$ 恆大於 0 ,使用位移操作 $d_m$ 到最後必為 0 ,因此可以 `d!=0` 作為中止條件。而 $c_{-1}=P_{0}2^{0}$ ,而 $P_0 = N$,因此 $c_{-1}$ 即為所求。
##### `ffs` / `fls` 取代 `__builtin_clz`
```c
int i_sqrt(int x)
{
if (x <= 1) /* Assume x is always positive */
return x;
int tmp = x, msb = 0;
while(tmp){
unsigned long i = __ffs(tmp);
msb += (i+1);
tmp >>= (i+1);
}
int z = 0;
for (int m = 1UL << (msb & ~1UL); m; m >>= 2) {
int b = z + m;
z >>= 1;
if (x >= b)
x -= b, z += m;
}
return z;
}
```
`__ffs` : find first bit in word
`__builtin_clz` : Returns the number of leading 0-bits in x
先使用一變數 `tmp` 等於 `x` ,由於 `__ffs` 在 lsb 一遇到非 0 bit 則停止,因此我使用一變數 msb 不斷累加 `__ffs` 回傳值 `i` ,同時`tmp` 不斷右移 `i` 直到 `tmp` 等於 0 ,不過須注意的是如果傳入 `__ffs` 之值最低 bit 為 1 時,將會回傳 0 而使 `tmp` 無法右移,將導致無窮迴圈,因此而外加上 1 ,代表會加上前方 bit 1 以避免此情況。
Linux 核心原始程式碼 [linux/lib/math/int_sqrt.c](https://github.com/torvalds/linux/blob/master/lib/math/int_sqrt.c)
```c
unsigned long int_sqrt(unsigned long x)
{
unsigned long b, m, y = 0;
if (x <= 1)
return x;
m = 1UL << (__fls(x) & ~1UL);
while (m != 0) {
b = y + m;
y >>= 1;
if (x >= b) {
x -= b;
y += m;
}
m >>= 2;
}
return y;
}
```
應用案例 [linux/block/blk-wbt.c](https://github.com/torvalds/linux/blob/486291a0e6246364936df1ecd64c90affef4b9c5/block/blk-wbt.c#L409)
```c
static void rwb_arm_timer(struct rq_wb *rwb)
{
struct rq_depth *rqd = &rwb->rq_depth;
if (rqd->scale_step > 0) {
/*
* We should speed this up, using some variant of a fast
* integer inverse square root calculation. Since we only do
* this for every window expiration, it's not a huge deal,
* though.
*/
rwb->cur_win_nsec = div_u64(rwb->win_nsec << 4,
int_sqrt((rqd->scale_step + 1) << 8));
} else {
/*
* For step < 0, we don't want to increase/decrease the
* window size.
*/
rwb->cur_win_nsec = rwb->win_nsec;
}
blk_stat_activate_nsecs(rwb->cb, rwb->cur_win_nsec);
}
```
這裡用途為透過 `request queue depth` 裡的成員 `scale_step` 動態決定 `cur_win_nsec` ,目的為調整執行速度及頻率,以便管理系統效能及資源。
調整方式為將 `win_nsec` 左移四位 再除上 `rqd->scale_step + 1) << 8` 之開跟號後的值。
#### 測驗 `2`
這裡要從最低位開始進行 mod 10 和 div 10 操作,不過因為一般整數格式不易取最低位值進行處理,因此以 char 格式作為輸入,再來看迴圈中
`int tmp = (b[i] - '0') + (a[i] - '0') + carry;` , ASCII 中,數字為 '0' 到 '9' ,扣除 '0' 範圍剩 '1' 到 '9' ,而 a 、 b 各一個再加上 carry ,因此最大為 19 ,也就是 `tmp` 不可能會大於的值。
除法因為 10 有包含 5 這個因數,無法完全用 2 的冪項來表示,因此結果會是不精確的。
因此僅能以接近 10 的倒數又同時符合直到小數點後第一位都是精確的條件下達成。
假設 $x$ 是除數,前面提到的 19 作範例,則:
$1.9\leq \frac{19}{x}\leq 1.99\\
\Rightarrow 9.55\leq x\leq10$
再來透過 bitwise operation 實作除法得到的除數 $q$ 必定是 $\frac{an}{2^N}$ 其中 $N$ 為非負整數, $a$ 正整數。因此只需要透過查看 $2^N$ 再配對適合的 $a$ 即可。
如:
$2^N=128, a=13,x = \frac{128}{13}\approx9.84$
因此除數 $q$ 可寫成 $\frac{13*n}{128}$
```c
unsigned d2, d1, d0, q, r;
d0 = q & 0b1;
d1 = q & 0b11;
d2 = q & 0b111;
q = ((((tmp >> 3) + (tmp >> 1) + tmp) << 3) + d0 + d1 + d2) >> 7;
r = tmp - (((q << 2) + q) << 1);
```
之所以為 `3`, `1`, `0` 是因為這些的 2 的冪可以組成 13 ($13 = 8 + 4 + 1 = 2^3 + 2^2 + 2^0$)
首先要用 bitwise operation 得到 $\frac{13tmp}{8} = \frac{tmp}{8}+\frac{tmp}{2}+tmp$
為得到 $13tmp$ ,因此左移 3 位。
`d0`, `d1`, `d2` 皆為保留右移後被捨棄的部分位元。
因此處選擇 $2^N=128$ ,將 `d0`, `d1`, `d2` 加上去後再右移 7 位做除以 128 的動作,得到除數 $q = \frac{13*n}{128}$。
餘數則透過 $tmp - (q \times 4 + q) \times 2 = tmp - 10\times{q}$ 。
包裝函式:
```c
*mod = in - ((q & ~0x7) + (*div << 1));
```
這裡 `q` 為右移三位 ($\frac{1}{8}$) 之前的 `div` ,又透過 `& ~0x7` 取第四位 ($2^3$) ,加上 `(*div << 1)` 得到 $10\times$*div。
##### 不依賴任何除法指令的 `% 9` (modulo 9) 和 `% 5` (modulo 5) 程式碼
> 實作不同於《Hacker's Delight》的程式碼
```c
void divmod_9(uint32_t in, uint32_t *div, uint32_t *mod)
{
while(in >= 9){
unsigned long high = __fls(in) - 3;
in -= (9 << high);
*mod += (1UL << high);
}
*div = in;
}
```
```c
void divmod_5(uint32_t in, uint32_t *div, uint32_t *mod)
{
while(in >= 5){
unsigned long high = __fls(in) - 2;
in -= (5 << high);
*mod += (1UL << high);
}
*div = in;
}
```
#### 測驗 `3`
ilog2 計算以 2 為底的對數
以二進位表示法在右移時計算同時是做 $/2$ 的動作,因此可以計算 log2 ,不過因為一次以一步的長度計算太慢,因此步數由 2 的冪開始。
不過這樣計算又使程式碼過於冗長,因此改以 `__builtin_clz` 直接計算頭部最後一個零在哪個位置以計算 log2 。
觀察 8 bits 二進位 `0b01001110` ,在頭部第一個 1 後不管有多少個 1 都不會造成此二進位數翻倍,因此不用計算後面有多少 1 。
Linux 核心原始程式碼 [linux/include/linux/log2.h](https://github.com/torvalds/linux/blob/master/include/linux/log2.h)
```c
int __ilog2_u32(u32 n)
{
return fls(n) - 1;
}
```
這裡有點類似上一題將 `__builtin_clz` 改成使用 `fls` 的 macro 。
應用案例 [linux/fs/ocfs2/ioctl.c](https://github.com/torvalds/linux/blob/712e14250dd2907346617eba275c46f53db8fae7/fs/ocfs2/ioctl.c#L402)
```c
static void o2ffg_update_histogram(struct ocfs2_info_free_chunk_list *hist,
unsigned int chunksize)
{
u32 index;
index = __ilog2_u32(chunksize);
if (index >= OCFS2_INFO_MAX_HIST)
index = OCFS2_INFO_MAX_HIST - 1;
hist->fc_chunks[index]++;
hist->fc_clusters[index] += chunksize;
}
```
這裡是利用 `__ilog2_u32` 尋找 `index` ,再對 `hist` 成員 `fc_chunks` 及 `fc_clusters` 的位置`index` 做特定處理。其為 [Oracle Cluster File System 2 (OCFS2)](https://en.wikipedia.org/wiki/OCFS2) (一般性用途的日誌檔案系統) 的直方圖更新函式。
> [ioctl](https://en.wikipedia.org/wiki/Ioctl) 系統呼叫的作用 參考 [Linux Ioctl internel](https://hackmd.io/@happy-kernel-learning/ByVFh6BEw)
當用 read/write 不能完成某一功能時,就用ioctl。
#### 測驗 `4`
觀察 EWMA 計算式
$S_t = \left\{
\begin{array}{l}
Y_0& t = 0 \\
\alpha Y_t + (1 - \alpha)\cdot S_{t-1}& t > 0 \\
\end{array}
\right.$
- $\alpha$ 表示歷史資料加權降低的程度,介在 0 ~ 1 之間,越高的 $\alpha$ 會使歷史資料減少的越快
- $Y_t$ 表示在時間 $t$ 時的資料點
- $S_t$ 表示在時間 $t$ 時計算出的 EWMA
$\alpha$ : `avg->weight`
$Y_t$ : `val`
$S_t$ : `avg->internal`
```c
struct ewma *ewma_add(struct ewma *avg, unsigned long val)
{
avg->internal = avg->internal
? (((avg->internal << avg->weight) - avg->internal) +
(val << avg->factor)) >> avg->weight
: (val << avg->factor);
return avg;
}
```
即
$S_t = \left\{
\begin{array}{l}
((\frac 1 \alpha - 1)\cdot S_{t-1} + Y_t)\alpha & t > 0 \\
Y_0& t = 0 \\
\end{array}
\right.$
就是上面給的數學式經過提取並變換順序。
Linux 核心原始程式碼 [linux/mm/ksm.c](https://github.com/torvalds/linux/blob/712e14250dd2907346617eba275c46f53db8fae7/mm/ksm.c#L378)
```c
/* Calculate exponential weighted moving average */
#define EWMA_WEIGHT 30
static unsigned long ewma(unsigned long prev, unsigned long curr)
{
return ((100 - EWMA_WEIGHT) * prev + EWMA_WEIGHT * curr) / 100;
}
```
Linux 核心原始程式碼 [linux/include/linux/average.h](https://github.com/torvalds/linux/blob/712e14250dd2907346617eba275c46f53db8fae7/include/linux/average.h#L9)
```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)); \
}
```
:::danger
分析程式碼的精確度,用實際數值帶入。
:::
可以發現基本上 `ewma` 與 `DECLARE_EWMA` 大同小異,差異僅為 `DECLARE_EWMA` 使用了位元運算,且固定 `weight` , `DECLARE_EWMA` 先將 factor 及 weight 取了 ilog2 以加速運算,而 `ewma` 因為應用於 KVM (花了很多時間在理解並寫下 [筆記](https://hackmd.io/@david96514/ryyZfLOeA#%E7%AC%AC-3-%E5%91%A8-%E6%B8%AC%E9%A9%97-4-%E6%89%80%E6%B5%AA%E8%B2%BB%E7%9A%84%E6%99%82%E9%96%93) ),掃描間隔太短,無法捨棄太多小數點,因此無法與 `DECLARE_EWMA` 使用相同寫法。
```c
/* Calculate exponential weighted moving average */
static unsigned long ewma(unsigned long prev, unsigned long curr)
{
return ((100 - EWMA_WEIGHT) * prev + EWMA_WEIGHT * curr) / 100;
}
```
以上兩邊都是 EWMA 不過計算方式明顯不同,<s>我也不太了解為什麼</s> 。
應用案例 [linux/drivers/net/wireless/mediatek/mt7601u/phy.c](https://github.com/torvalds/linux/blob/master/drivers/net/wireless/mediatek/mt7601u/phy.c#L970)
```c
static void mt7601u_agc_tune(struct mt7601u_dev *dev)
{
u8 val = mt7601u_agc_default(dev);
long avg_rssi;
if (test_bit(MT7601U_STATE_SCANNING, &dev->state))
return;
/* Note: only in STA mode and not dozing; perhaps do this only if
* there is enough rssi updates since last run?
* Rssi updates are only on beacons and U2M so should work...
*/
spin_lock_bh(&dev->con_mon_lock);
avg_rssi = ewma_rssi_read(&dev->avg_rssi);
spin_unlock_bh(&dev->con_mon_lock);
if (avg_rssi == 0)
return;
avg_rssi = -avg_rssi;
if (avg_rssi <= -70)
val -= 0x20;
else if (avg_rssi <= -60)
val -= 0x10;
if (val != mt7601u_bbp_rr(dev, 66))
mt7601u_bbp_wr(dev, 66, val);
/* TODO: also if lost a lot of beacons try resetting
* (see RTMPSetAGCInitValue() call in mlme.c).
*/
}
```
參考 [IoT 的基礎知識 - RSSI](https://peterpowerfullife.com/blog/tech-rssi/)
> RSSI 是衡量設備從接收端接收信號能力的指標
參考 [UniFi Network - 了解和使用最小 RSSI](https://help.tw.ui.com/articles/221321728/)
> 主要目的是幫助客戶端在兩個臨近的 AP (Assess Point) 之間漫遊。它可以防止設備在較弱的信號強度下“卡住”,一直連接初始 AP,而不會漫遊到性能更佳的新 AP。一旦信號低於設置的最小 RSSI 值,初始 AP 會踢掉客戶端,以便它可以重新連接到新 AP。
而 EWMA 使經過時間越久的歷史資料的權重也會越低,用 EWMA 將過去一段時間信號強度乘上不同的權值來計算評估信號強度,若發現信號強度不足,則變更 AP 的選取。
#### 測驗 `5`
先透過 x 得到最高位 bit 再透過遞減逐漸逼近正確 log2 值。
`(r | shift | x > 0x2) + 1`
`r` 負責 bit >= 2 , `shift` 負責 bit = 1 , `x > 0x2` 則負責 bit = 0。
最後加 1 則是因為 ceil 。
##### 改進程式碼,使其得以處理 `x = 0` 的狀況,並仍是 branchless
```c
int ceil_ilog2(uint32_t x)
{
uint32_t r, shift;
int mask = (x == 0) << 3;
x--;
r = (x > 0xFFFF) << 4;
x >>= r;
shift = (x > 0xFF) << 3;
x >>= shift;
r |= shift;
shift = (x > 0xF) << 2;
x >>= shift;
r |= shift;
shift = (x > 0x3) << 1;
x >>= shift;
return ((r | shift | x > 1) + 1) >> mask;
}
```
這裡我想法是若 `x==0` 則將 1 右移 3 位 (只要數值大於 2 就好,因為原本輸入為 0 時輸出結果為 32 ,因此需右移 5 位,而 1 右移兩位為 4 ,不足以將 32 歸 0) ,最後透過右移 mask 將 32 歸 0 (原不等於 0 的值不會受影響)。
Linux 核心原始程式碼 [linux/tools/include/linux/log2.h](https://github.com/torvalds/linux/blob/master/tools/include/linux/log2.h#L151)
```c
#define roundup_pow_of_two(n) \
( \
__builtin_constant_p(n) ? ( \
(n == 1) ? 1 : \
(1UL << (ilog2((n) - 1) + 1)) \
) : \
__roundup_pow_of_two(n) \
)
```
應用案例 [linux/tools/perf/util/print_binary.c](https://github.com/torvalds/linux/blob/39cd87c4eb2b893354f3b850f916353f2658ae6f/tools/perf/util/print_binary.c#L15)
```c
int binary__fprintf(unsigned char *data, size_t len,
size_t bytes_per_line, binary__fprintf_t printer,
void *extra, FILE *fp)
{
size_t i, j, mask;
int printed = 0;
if (!printer)
return 0;
bytes_per_line = roundup_pow_of_two(bytes_per_line);
mask = bytes_per_line - 1;
printed += printer(BINARY_PRINT_DATA_BEGIN, 0, extra, fp);
for (i = 0; i < len; i++) {
if ((i & mask) == 0) {
printed += printer(BINARY_PRINT_LINE_BEGIN, -1, extra, fp);
printed += printer(BINARY_PRINT_ADDR, i, extra, fp);
}
printed += printer(BINARY_PRINT_NUM_DATA, data[i], extra, fp);
if (((i & mask) == mask) || i == len - 1) {
for (j = 0; j < mask-(i & mask); j++)
printed += printer(BINARY_PRINT_NUM_PAD, -1, extra, fp);
printer(BINARY_PRINT_SEP, i, extra, fp);
for (j = i & ~mask; j <= i; j++)
printed += printer(BINARY_PRINT_CHAR_DATA, data[j], extra, fp);
for (j = 0; j < mask-(i & mask); j++)
printed += printer(BINARY_PRINT_CHAR_PAD, i, extra, fp);
printed += printer(BINARY_PRINT_LINE_END, -1, extra, fp);
}
}
printed += printer(BINARY_PRINT_DATA_END, -1, extra, fp);
return printed;
}
```
:::danger
file 的翻譯是「檔案」,而非「文件」(document)
:::
目的是以指定的格式將二進位資料輸出到檔案中。
這段函式雖然用到 ceil 和 log2 的組合,不過我不知道為什麼要這樣做,唯一的想法是為了要控制在 2 的冪以利用位元運算達成更快的計算是否達到該行結尾。
### 第 4 周測驗題
#### 測驗 `1`
[漢明距離](https://en.wikipedia.org/wiki/Hamming_distance)
即兩個等長字串之間的 Hamming distance ,是兩個字串對應位置的不同字符的個數。
而計算方法為 `__builtin_popcount(x ^ y)`
```c
int totalHammingDistance(int* nums, int numsSize)
{
int total = 0;
for (int i = 0;i < numsSize;i++)
for (int j = 0; j < numsSize;j++)
total += __builtin_popcount(nums[i] ^ nums[j]);
return total >> 1;
}
```
若輸入有 n 個字串,則這 n 個字串會與所有 n 個字串做 Hamming distance 的計算。因此任意兩的字串會以不同的順序各被計算一次,因會被重複計算,要右移一位以得到正確值。
寫出避免重複運算的版本。
```c
int totalHammingDistance(int* nums, int numsSize)
{
int total = 0;
for (int i = 0;i < numsSize;i++)
for (int j = i+1; j < numsSize;j++)
total += __builtin_popcount(nums[i] ^ nums[j]);
return total;
}
```
##### 不依賴 popcount,針對 LeetCode 477. Total Hamming Distance 撰寫出更高效的程式碼。
> 所有 Number 中從第 0 個 bit 開始計算全部該 bit 的 Hamming Distance
因此目標被設定為時間複雜度 $O(n)$ 因 bit 數不被輸入數量影響,只與作業系統有關。
```c
int totalHammingDistance(int* nums, int numsSize)
{
int total = 0;
int num_1s = 0;
int displace = 0;
int bits = 32;
for (int i = 0;i < bits;i++){
for (int j = 0; j < numsSize;j++){
num_1s += ((nums[j]>>displace) & 0x1);
-(~displace);
}
total += num_1s * (numsSize - num_1s);
}
return total;
}
```
#### 測驗 `2`
$2^k \equiv \begin{cases}
1 (mod \ \ 3), \ \ k \ \ even\\
-1 (mod \ \ 3), \ \ k \ \ odd\\
\end{cases}$
n 的 2 進位表示法 $b_{n-1}b_{n-2}\cdots b_1 b_0$
$n \equiv \sum_{i=0}^{n-1} b_i\cdot (-1)^i \ \ (mod \ \ 3)$
透過 `n = popcount(n & 0x55555555) - popcount(n & 0xAAAAAAAA)` 計算位元和, 5 的二進位表示為 0101 (偶數位置 0 、 2 ),而 A 為 1010 (奇數位置 1 、 3 )
再透過定理
$popcount(x \And \overline{m}) - popcount(x \And m) = popcount(x \bigoplus m) - popcount(m)$
將
`n = popcount(n & 0x55555555) - popcount(n & 0xAAAAAAAA)`
變為
`n = popcount(n ^ 0xAAAAAAAA) - 16`
```c
int mod3(unsigned n)
{
n = popcount(n ^ 0xAAAAAAAA) + 23;
n = popcount(n ^ 0x2A) - 3;
return n + ((n >> 31) & 3);
}
```
`+ 23`:原先是 `-16` ,這裡因為避免 n 變為負數,因此加上一足夠大的 3 的倍數。
接著重複應用於得到的 n 上直到 n 介於 0~2 (類似不斷除以 3)。
最後為了避免 n 為負,將 n 右移 31 位(原本的正數則不受影響)並 `& 3` 再加上原本的 n 。
**tic-tac-toe**
#### 測驗 `3`