# 2020q3 Homework (quiz5)
contributed by < `shauming1020` >
###### tags: `homework`
## 測驗 1
```c=
#include <stdio.h>
#include <stdlib.h>
double divop(double orig, int slots) {
if (slots == 1 || orig == 0)
return orig;
int od = slots & 1;
double result = divop(orig / D1, od ? (slots + D2) >> 1 : slots >> 1);
if (od)
result += divop(result, slots);
return result;
}
```
> 假設 divop() 的第二個參數必為大於 0 的整數，而且不超過 int 型態能表達的數值上界。
### １. 解釋上述程式碼運作原理
- ```int od = slots & 1;``` 判斷 slots 是否為奇數
- ```double result = divop(orig / D1, od ? (slots + D2) >> 1 : slots >> 1);```
- ```od ? (slots + D2) >> 1 : slots >> 1)``` 當 slots 是奇數時，必須 + 1 後再進行右移 (/2) ，如此才會滿足遞迴中止條件 (slots == 1)，因此 D2 為 1
- ```orig / D1``` 浮點數無法透過右移 bit 來進行除 2 的動作，因此 D1 為 2
- 將 orig 表示為分數，分子與分母同除 2 後不會影響最後結果，當分母為 1 時，分子即為商數
- ```if (od) result += divop(result, slots);```
- 當 slots 為奇數時，slots 加 1 後再去除 2 會產生誤差，導致求出的商數比原商數小
> e.g. $1 \div 3 = 0.333$ 在計算時會變成 $1 \div 4 = 0.250$，誤差為 $0.083$，恰近似於 $0.250 / 3$
> 將上述例子用數學式描述 $\dfrac {1}{3} = \dfrac{1}{4} + \dfrac{\dfrac{1}{4}}{3}=\dfrac{1}{4}+\dfrac{1}{12}=\dfrac{4}{12}$
- 因此需要將誤差給加上去，而 slots 為偶數時，可以直接同除 2，不會有誤差的產生
- 所以每次求誤差時會一直遞迴執行程式，e.g. 求 $\dfrac{\dfrac{1}{4}}{3}$ 時會再去計算誤差
### 2. 以編譯器最佳化的角度，推測上述程式碼是否可透過 tail call optimization 進行最佳化，搭配對應的實驗;
> 搭配閱讀 C 語言: 編譯器和最佳化原理 及 C 語言: 遞迴呼叫
觀察 gcc -S -O 編譯器行為
- Recursive
```c=
int fib(int n) {
if (n == 0) return 0;
if (n == 1) return 1;
return fib(n - 1) + fib (n - 2);
}
main()
{
fib(100);
}
```
- Assembly code
```ass
...
fib:
.LFB0:
.cfi_startproc
pushq %rbp # Push to stack
.cfi_def_cfa_offset 16
.cfi_offset 6, -16
pushq %rbx # Push to stack
.cfi_def_cfa_offset 24
.cfi_offset 3, -24
subq $8, %rsp
.cfi_def_cfa_offset 32
movl %edi, %ebx
testl %edi, %edi
je .L2
cmpl $1, %edi
jne .L4
.L2:
movl %ebx, %eax
addq $8, %rsp
.cfi_remember_state
.cfi_def_cfa_offset 24
popq %rbx
.cfi_def_cfa_offset 16
popq %rbp
.cfi_def_cfa_offset 8
ret
.L4:
.cfi_restore_state
leal -1(%rdi), %edi
call fib
movl %eax, %ebp
leal -2(%rbx), %edi
call fib
leal 0(%rbp,%rax), %ebx
jmp .L2
.cfi_endproc
...
main:
.LFB1:
.cfi_startproc
subq $8, %rsp
.cfi_def_cfa_offset 16
movl $100, %edi
call fib
movl $0, %eax
addq $8, %rsp
.cfi_def_cfa_offset 8
ret
.cfi_endproc
...
```
- Tail recursion
```c=
int fib(int n, int a, int b) {
if (n == 0) return a;
return fib(n - 1 , b, a + b);
}
main()
{
fib(100, 0, 1);
}
```
---
## 測驗 2
假設 float 為 32-bit 寬度，考慮以下符合 IEEE 754 規範的平方根程式，請補上對應數值，使其得以運作：
```c=
#include <stdint.h>
/* A union allowing us to convert between a float and a 32-bit integers.*/
typedef union {
float value;
uint32_t v_int;
} ieee_float_shape_type;
/* Set a float from a 32 bit int. */
#define INSERT_WORDS(d, ix0) \
do { \
ieee_float_shape_type iw_u; \
iw_u.v_int = (ix0); \
(d) = iw_u.value; \
} while (0)
/* Get a 32 bit int from a float. */
#define EXTRACT_WORDS(ix0, d) \
do { \
ieee_float_shape_type ew_u; \
ew_u.value = (d); \
(ix0) = ew_u.v_int; \
} while (0)
static const float one = 1.0, tiny = 1.0e-30;
float ieee754_sqrt(float x)
{
float z;
int32_t sign = 0x80000000;
uint32_t r;
int32_t ix0, s0, q, m, t, i;
EXTRACT_WORDS(ix0, x);
/* take care of zero */
if (ix0 <= 0) {
if ((ix0 & (~sign)) == 0)
return x; /* sqrt(+-0) = +-0 */
if (ix0 < 0)
return (x - x) / (x - x); /* sqrt(-ve) = sNaN */
}
/* take care of +INF and NaN */
if ((ix0 & 0x7f800000) == 0x7f800000) {
/* sqrt(NaN) = NaN, sqrt(+INF) = +INF,*/
return x;
}
/* normalize x */
m = (ix0 >> 23);
if (m == 0) { /* subnormal x */
for (i = 0; (ix0 & 0x00800000) == 0; i++)
ix0 <<= 1;
m -= i - 1;
}
m -= 127; /* unbias exponent */
ix0 = (ix0 & 0x007fffff) | 0x00800000;
if (m & 1) { /* odd m, double x to make it even */
ix0 += ix0;
}
m >>= 1; /* m = [m/2] */
/* generate sqrt(x) bit by bit */
ix0 += ix0;
q = s0 = 0; /* [q] = sqrt(x) */
r = QQ0; /* r = moving bit from right to left */
while (r != 0) {
t = s0 + r;
if (t <= ix0) {
s0 = t + r;
ix0 -= t;
q += r;
}
ix0 += ix0;
r >>= 1;
}
/* use floating add to find out rounding direction */
if (ix0 != 0) {
z = one - tiny; /* trigger inexact flag */
if (z >= one) {
z = one + tiny;
if (z > one)
q += 2;
else
q += (q & 1);
}
}
ix0 = (q >> 1) + QQ1;
ix0 += (m << QQ2);
INSERT_WORDS(z, ix0);
return z;
}
```
### 1. 解釋上述程式碼何以運作
> 搭配研讀 [以牛頓法求整數開二次方根的近似值](http://science.hsjh.chc.edu.tw/upload_works/106/44f6ce7960aee6ef868ae3896ced46eb.pdf) (第 19 頁)
---
## 測驗 3
$kx = N - \dfrac{(k-1)k}{2}$，由於 x 要為正整數，因此 $x = (N - \dfrac{(k-1)k}{2}) / k \in Z^+$，即 $(N - \dfrac{(k-1)k}{2}) \% k = 0$，而 k 的範圍可從近似解 $k < \sqrt{2N}$ 得知
```c=
int consecutiveNumbersSum(int N)
{
if (N < 1)
return 0;
int ret = 1;
int x = N;
for (int i = 2; i < x; i++) {
x += ZZZ;
if (x % i == 0)
ret++;
}
return ret;
}
```
### 1. 解釋上述程式碼何以運作
- 撰寫與推導後的數學式相同的程式碼
```c=
int consecutiveNumbersSum(int N)
{
if (N < 1)
return 0;
int ret = 1;
int x = N;
for (int i = 2; i < sqrt(2*N); i++) {
if ((N - (i-1) * i / 2) % i == 0) ret++;
}
return ret;
}
```
- 觀察 $S_k = N - \dfrac{(k-1)k}{2}, 1 \le k < \sqrt{2N}$
|$k$| $S_k$ | $S_k - S_{k-1}$ |
| --| -------- | --------------- |
| 1 | N | 0 |
| 2 | N + (-1) | -1 |
| 3 | N + (-3) | -2 |
| 4 | N + (-6) | -3 |
| 5 | N + (-10)| -4 |
| n | N + $\dfrac{(-n+1)n}{2}$| 1-n |
- 檢驗 $S_k \% k = 0$ 是否成立來判斷正整數解 x，從上表可以觀察出 $S_k = S_{k-1} + (1-k)$，故 ZZZ 為 (1 - i)
> e.g. $S_5 = N + (-10) = S_4 + (1 - 5) = N +(-6) + (-4)$
- 而 $S_k \% k = 0$，且 k 為大於等於 2 的正整數，因此 $2 \le k \le S_k$
- 在 for 迴圈內的 ```x += (1-i) ``` 計算出 $S_k$ 後，**進入下個迴圈前，k 會先被加一再和 $S_k$ 比較大小**，因此範圍為 $(k+1) < S_k$，迴圈的條件判斷式才會是 ``` i < x ```
### 2. 嘗試改寫為更有效率的實作
---
Sign in with Wallet
Connect another wallet