# 2020q3 Homework5 (quiz5)
contributed by < `Holycung` >
[2020q3 Homework5 (quiz5) 題目](https://hackmd.io/@sysprog/2020-quiz5)

## 目錄
[TOC]

## 測驗 1
### 題目
考慮到以下浮點數除法程式: (fdiv.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 型態能表達的數值上界。請補完程式碼。

### 作答
回答這題要先了解一下在 [IEEE 754](https://zh.wikipedia.org/wiki/IEEE_754) 下浮點數的定義。


這題的觀念是把原本的浮點數除法都替換成除以 2 進行約分，因為整數除以二可以用往右 shift 一位來加速，浮點數除以二只需要把 exponent 的部分減一就好。那就可以把狀況分成下面兩種，假設 $X \div Y$，已知 Y 為大於零的整數，第一種狀況就是 $Y$ 為偶數，那我們可以直接對 $X\ Y$ 都向右 shift 一位進行約分，也就是 $X \div Y = \frac{X}{2} \div \frac{Y}{2}$。第二種狀況就是 $Y$ 為奇數，這就是我們所要探討的。

$Y$ 為奇數了話我們就沒有辦法直接除二，那這邊的做法就是先將 $Y+1$，把原本的 $X \div Y$ 變成 $X \div (Y+1)$，這樣 $Y+1$ 就是偶數就可以用第一種狀況的方法，不過我們知道 $X \div (Y+1)$ 跟原本的值會不同，且 $X \div Y > X \div (Y+1)$，所以我們會需要補差值回去，那個**差值**就是 $X \div Y - X \div (Y+1)$，稍微化減一下如下。

\begin{aligned}
差值 &= \frac{X}{Y} - \frac{X}{Y+1} \\
&= \frac{X(Y+1)}{Y(Y+1)} - \frac{XY}{Y(Y+1)} \\
&= \frac{X}{Y(Y+1)}
\end{aligned}

所以這邊我們得到了**差值** $=\frac{X}{Y(Y+1)}$，可以發現這邊 $\frac{X}{(Y+1)}$ 跟我們前面第二種狀況是一樣的，接下來我們就可以來看程式。

#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;
}

第 5 行這邊是判斷如果，除數 `slot` 是 1 或是被除數 `orig` 是 0 就直接回傳被除數 `orig`。
第 7 行 `od` 是判斷除數是否為奇數，只要拿最後一個位元跟 1 做 bitwise and 就可以得到。
第 8 行這邊就是我們剛剛上面推論的方法，把浮點數的除法全部換成除以二約分遞迴下去，`divop` 第一個參數是 `orig / D1`，就是對被除數除以二的動作，所以 `D1` = `2`這邊是浮點數除二所以能用 shift。
第二個參數會先判斷除數 `slots` 是否為奇數，如果是偶數就直接進行 shift，如果是奇數就加一後在進行 shift，所以 `D2` = `1`。最後把遞迴的結果存回到 `result` 變數，所以 `result` = $\frac{X}{(Y+1)}$。
第 9 10 行最後就是要判斷除數是否為奇數，如果是奇數了話我們需要把**差值**補回去，經過剛剛的推導**差值** $=\frac{X}{Y(Y+1)}$，也就是 $\frac{result}{Y}$。

### 延伸問題
:::info
TODO
- 以編譯器最佳化的角度，推測上述程式碼是否可透過 tail call optimization 進行最佳化，搭配對應的實驗;
- 比照 浮點數運算和定點數操作，改寫上述程式碼，提出效率更高的實作，同樣避免使用到 FPU 指令 (可用 gcc -S 觀察對應的組合語言輸出)，並與使用到 FPU 的實作進行效能比較
:::

## 測驗 2
### 題目
延續 從 [$√2$ 的運算談浮點數](https://hackmd.io/@sysprog/ryOLwrVtz?type=view)，假設 float 為 32-bit 寬度，考慮以下符合 IEEE 754 規範的平方根程式，請補上對應數值，使其得以運作：

#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;
}

請補完程式碼，使上述程式碼得以運作。

### 作答
再回答這題前可以先看過[浮點數運算和定點數操作](https://hackmd.io/@NwaynhhKTK6joWHmoUxc7Q/H1SVhbETQ)，熟悉一下 IEEE 754 對浮點數的規範，再來直接看到程式碼。

/* 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;

在 C 當中的 union 用法可以參考 [Union type wiki](https://en.wikipedia.org/wiki/Union_type)。
> a union is a value that may have any of several representations or formats within the same position in memory;

/* 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)

這邊有用到 [你所不知道的C語言：技巧篇](https://hackmd.io/@sysprog/c-trick?type=view#do----while0-%E5%B7%A8%E9%9B%86) 提到的 `do { ... } while(0) 巨集` 是為了避免 dangling else 問題，也可以參考這篇 [macro do{}while(0)](http://github.tiankonguse.com/blog/2014/09/30/do-while.html) 更詳細的解釋。
這兩個 macro `INSERT_WORDS`、`EXTRACT_WORDS`，是運用 union 來轉換 float 跟 int，參考 [RinHizakura](https://hackmd.io/@RinHizakura/SJnXiSewD) 的例子如下:

SIGN  EXPONENT  FRACTION
0     00000010  01000000000000000000000

如果用 `uint32_t` 的角度來解釋，就是 $2^{24}+2^{21}$
如果從 `float` (IEEE 754) 的角度來解釋，則是 $(1+1*2^{-2})×2^{(2−127)}$

    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 */
    }

第 34 行先把浮點數 x 轉成整數 `ix0`。
第 37 行這邊是先處理小魚等於 0 的數字，第 38 行是判斷是否為正零或負零(除了 sign bit 其他都是 0)，如果是 0 的話就回傳 0。
第 40 行是判斷是否小於 0，對小於零的數取平方根的話要回傳 NaN，第 41 行是透過 0.0 / 0.0 產生 NaN，這邊可以從 [IEEE 754-2019](https://ieeexplore.ieee.org/stamp/stamp.jsp?tp=&arnumber=8766229) 規格書當中 6.2.1 先了解到 NaN 的編碼格式，然後又分程 quiet NaN 和 signal NaN，再看到 7.2 有明確規範當浮點數 0.0 / 0.0 的時候會產生 quiet NaN。
> **6.2.1 NaN encodings in binary interchange formats**
> All binary NaN bit strings have the sign bit S set to 0 or 1 and all the bits of the biased exponent field E set to 1 (see 3.4). A quiet NaN bit string should be encoded with the first bit (d1) of the trailing significand field T being 1. A signaling NaN bit string should be encoded with the first bit of the trailing significand field being 0.

> **7.2 Invalid operation**
> For operations producing results in floating-point format, the default result of an operation that signals the invalid operation exception shall be a quiet NaN that should provide some diagnostic information.
> e) division: division(0, 0) or division(∞, ∞)

    /* take care of +INF and NaN */
    if ((ix0 & 0x7f800000) == 0x7f800000) {
        /* sqrt(NaN) = NaN, sqrt(+INF) = +INF,*/
        return x;
    }

第 44 行是判斷 +INF 和 NaN 如果是，`ix0 & 0x7f800000 == 0x7f800000` 就可以判斷 Exponent 是否等於 255，也就是 INF、NaN 的狀況。


    /* 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] */

第 48 行把 ix0 向右 shift 23 個得到 m 就是 exponent 的值，如果 `m == 0` 代表是 denormalized number，這邊就需要先把他變成 normalized number 的形式一起處理，所以這邊的想法是把 ix0 往左 shift 也就是把 fraction 的部分往左移動，直到第 24 個位元是 1，然後記錄下移動了幾次 `i`，每移動一次就相當指數 m 要減一，所以最後要把 `m -= i - 1`，就變成了 normalized 的形式。

因為剛剛處理完 denormalized，這邊就可以直接用 normalized 的形式處理，也就是 $1.F \times 2^{E-127}$，所以第 55 行就把 m 減掉 127。

第 56 行是把 `ix0` 只留下後面 fraction 23 bits 的部分，把 24 位元變成 1，因為前面已經把 exponent 的部分保留下來到 m，然後也把 denormalized 的部分變成 normalized 的形式，所以第 24 位元本來就是 1，因此這邊就只留下後面 24 位元也就是 1.F。

第 57 行是判斷 m 是否為偶數，因為對 $2^E$ 取平方根了話就會是 $2^{E/2}$，那如果 m 不是偶數了話，就把 `ix0` 乘以二，這樣就可以把 m 除以二。

    /* 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;      // t = s_i + r_i
        if (t <= ix0) {      
            s0 = t + r;  // s_{i+1} = t + r_i
            ix0 -= t;    // x_{i+1}/2 = x_i - t
            q += r;      // q_{i+1} = q_i + r_i
        }
        ix0 += ix0;      // x_{i+1}
        r >>= 1;
    }

接著為計算 square root 的部分，這邊的方法可以參考 [How to calculate a square root without a calculator]( https://www.homeschoolmath.net/teaching/square-root-algorithm.php)，其方法的證明可以參考 [Why the square root algorithm works](https://www.homeschoolmath.net/teaching/sqr-algorithm-why-works.php)，應用到二進位 [Methods of computing square roots](https://en.wikipedia.org/wiki/Methods_of_computing_square_roots#Binary_numeral_system_(base_2)) 也是一樣的，我們則可以利用這個方法把平方根換成較快的 bit shift operations，接下來進行推導，推導這部分是參考 [RusselCK](https://hackmd.io/@RusselCK/sysprog2020_quiz5#Q2-%E6%B5%AE%E9%BB%9E%E6%95%B8%E9%96%8B%E4%BA%8C%E6%AC%A1%E6%96%B9%E6%A0%B9%E7%9A%84%E8%BF%91%E4%BC%BC%E5%80%BC)。

- 這邊我們要計算 $x$ 的平方根 $\sqrt{x}$ bit by bit 的數學推導。
- 令 $q_i=\sqrt{x},\ q_i$ 代表到小數點後面 $i$ 位精確值 $i \ge 0$
- 那假設今天 $1 \le x < 4$，那 $\sqrt{x}$ 的整數位必定是 1，所以 $q_0=1$，就會後面會探討為甚麼要這樣假設。
- 考慮 $q_{i+1}$，也就是 bit by bit 的逼近下去:
1. 若 $(q_{i} + 2^{-(i+1)})^2 \le x$，則 $q_{i+1}=q_i+2^{-(i+1)}$
- 即 $q_i$ 下一位補 1
2. 若 $(q_{i} + 2^{-(i+1)})^2 > x$，則 $q_{i+1}=q_i$
- 即 $q_i$ 下一位補 0
- 整理上面第一種情況 $(q_{i} + 2^{-(i+1)})^2 \le x$：
$$
\begin{aligned}
{(q_i + 2^{-(i+1)})}^2 \leq x &\Rightarrow q_i^2 + 2\times q_i \times (2^{-(i+1)}) + (2^{-(i+1)})^2 \leq x \\
&\Rightarrow q_i \times (2^{-i}) + 2^{-2(i+1)} \leq x - q_i^2 \\
&\Rightarrow q_i \times 2 \space + 2^{-(i+1)} \leq (x - q_i^2) \times 2^{(i+1)} \\
&\Rightarrow s_i + r_i \leq x_i, \\
&\begin{split} where \space &s_i = 2 \times q_i, \\
&r_i = 2^{-(i+1)}, \\
&x_i = (x - q_i^2) \times 2^{(i+1)} = (x - q_i^2) \times r_i^{-1}
\end{split}
\end{aligned}
$$
- 因此
- 若 $s_i + r_i \leq x_i$，則
* $q_{i+1}= q_i + r_i$
* $s_{i+1}=(s_i + r_i) + r_i$
* $x_{i+1}=2x_i - 2(s_i + r_i) = 2[x_i - (s_i + r_i)]$
:::spoiler 推導過程
$$
\begin{split}
s_{i+1} &= 2 \times q_{i+1} \\
&= 2 \times [q_i + 2^{-(i+1)}] \\
&= (2 \times q_i) + 2^{-i} \\
&= s_i + 2^{-i} \\
&= s_i + 2^{-(i+1)} + 2^{-(i+1)} \\
&= (s_i + r_i) + r_i \\
\end{split}
$$
$$
\begin{split}
x_{i+1} &=[x - q_{i+1}^2] \times 2^{((i+1)+1)} \\
&= [x - (q_i + r_i)^2] \times 2^{(i+2)} \\
&= [x - (q_i^2 + 2q_ir_i + r_i^2)] \times 2r_i^{-1} \\
&= [x - q_i^2] \times 2r_i^{-1} - 2q_ir_i \times 2r_i^{-1} - r_i^2 \times 2r_i^{-1} \\
&= 2x_i - 2s_i - 2r_i \\
&= 2x_i - 2(s_i + r_i) \\
\end{split}
$$
:::
- 若 $s_i + r_i > x_i$，則
* $q_{i+1} = q_i$
* $s_{i+1}=s_i$
* $x_{i+1}=2x_i$
:::spoiler 推導過程
$$
\begin{split}
s_{i+1} &= 2 \times q_{i+1} \\
&= 2 \times q_i \\
&= s_i \\
\end{split}
$$
$$
\begin{split}
x_{i+1} &=[x - q_{i+1}^2] \times 2^{((i+1)+1)} \\
&= [x - q_i^2] \times 2^{(i+2)} \\
&= [x - q_i^2] \times 2r_i^{-1} \\
&= 2x_i
\end{split}
$$
:::

接下來我們可以回到程式碼。

    /* 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;      // t = s_i + r_i
        if (t <= ix0) {      
            s0 = t + r;  // s_{i+1} = t + r_i
            ix0 -= t;    // x_{i+1}/2 = x_i - t
            q += r;      // q_{i+1} = q_i + r_i
        }
        ix0 += ix0;      // x_{i+1}
        r >>= 1;
    }

第 63 行，把 `ix0` 往左 shift 一位。這樣 `ix0` 的寬度可能是 25 或 26 個位元，所以 `QQ0` = `0x01000000`，第 25 個位元為一，等於把小數點點在第 24 25 個位元之間，從第 25 個 bit 開始，這樣代表我們會精確到小數點後面第 24 位，那在 IEEE 754 的規範中 fraction 部分有 23 位，我們就可以把最後一位做 rounding，精確到小數點後 23 位。小數點點在這邊也可以確保我們前面的假設 $1 \le x < 4$ 是正確的。

static const float one = 1.0, tiny = 1.0e-30;

    /* 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);
        }
    }

第 79 行這邊是判斷 rounding，首先要先判斷系統的 rounding mode，這邊透過 1 加減一個很小的數 `tiny`，讓系統進行 rounding 來判斷其 rounding mode。
- 先用 `one - tiny` 如果小於 1，代表是 rounding down，那不管最低位是 0 或 1 都會被捨棄掉，所以不需要額外動作。
- 接下來判斷 `one + tiny` 如果大於 1，代表是 rounding up，那 23 bit 就直接直接進位，所以這邊是 `q += 2`。
- 最後如果 `one + tiny` 如果等於 1，代表是 round to nearest，此時就要判斷最低位是 0 或 1，如果是 1 就要進位，所以是 `q += (q & 1)`。

    ix0 = (q >> 1) + QQ1;
    ix0 += (m << QQ2);
	
    INSERT_WORDS(z, ix0);

    return z;

最後要把 exponent 跟 fraction 的部分合併回來，所以第 90 行先把 q 往右 shift 1 位，不過要注意 q 現在還剩下 24 位，因為剛剛整理得到的是 $1.F$ 所以我們這邊的 `QQ1` = `0x3f000000`，讓中間 exponent 最低 7 位 bit 都補上一，因為 IEEE 754 normalized 的格式是 $1.F \times 2^{E-127}$，也就是在這邊 exponent 先加上 127。第 91 行，在把 m 往左 shift 23 個，加回去 `ix0`。

最後在透過 `INSERT_WORDS` 轉回浮點數回傳。

### 延伸問題
練習 LeetCode [69. Sqrt(x)](https://leetcode.com/problems/sqrtx/)，提交不需要 FPU 指令的實作。

#### 作答
[從 $√2$ 的運算談浮點數](https://hackmd.io/@sysprog/ryOLwrVtz?type=view) 了解到開根號可以透過二分逼近法或是十分逼近法來找到，再研讀 [以牛頓法求整數開二次方根的近似值 (第 19 頁)](http://science.hsjh.chc.edu.tw/upload_works/106/44f6ce7960aee6ef868ae3896ced46eb.pdf)，決定使用牛頓法來實作。

$\sqrt[n]{N} ≒ [(n-1) \times a + \frac{N}{a^{n-1}}] \div n = a_1$

透過多次逼近可以取得更精準的值，今天我們要找到平方根，所以 n=2:

$\sqrt{N} ≒ [(2-1) \times a + \frac{N}{a^{2-1}}] \div 2 = \dfrac{a + \frac{N}{a}}{2}$

接下來看到程式。

int mySqrt(int x){
    if (x == 0 || x == 1)
        return x;
    if (x == INT_MAX) // avoid overflow
        x-=1;
    
    unsigned int a; // using right shift
    a = (x + 1) >> 1;
    
    while (a > x/a){
        a = (a + x/a) >> 1;
    }
    return a;
}

那最開始的值要從多少開始逼近呢，這邊可以從 INT_MAX 最大值開始逼近，但是我們可以透過算幾不等式，找到比 $\sqrt x$ 大但又不會太大的值開始，轉換成程式碼就是 `(x + 1) >> 1`，不過這邊要小心因為會對 x 做加一的動作，所以如果 `x == INT_MAX` 要避免 overflow 就要先把他減一。

\begin{align}
\dfrac{a+b}{2} \ge& \sqrt{ab},\ where\ a,\ b>0 \\
=> \dfrac{a+1}{2} \ge& \sqrt{a},\ where\ a>0,\ b=1
\end{align}

最後 while 迴圈的終止條件則是找到第一個 $a^2 \le x$ 的 $a$。


---
另外還有找到一個[快速整數平方根算法](http://www.azillionmonkeys.com/qed/ulerysqroot.pdf)，有時間了話再來研究QQ

int mySqrt(int v)
{
    unsigned long temp, nHat = 0, b = 0x8000, bshft = 15;
    do{
        if (v >= (temp = (((nHat << 1) + b) << bshft--)))
        {
            nHat += b;
            v -= temp;
        }
    } while (b >>= 1);
    return nHat;
}

## 測驗 3
### 題目
LeetCode [829. Consecutive Numbers Sum](https://leetcode.com/problems/consecutive-numbers-sum/) 給定一個正整數 N，問 N 能寫成多少種連續正整數之和，例如 9 可寫成 4+5，或者　2+3+4。由於要寫成連續正整數之和，則可用等差數列來思考，且差值為 1，這個等差數列不必從 1 開始。
參考實作如下:

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;
} 

請補完程式碼，使上述程式碼得以運作。

### 作答
假設從 x 開始，且個數共有 k 個，則可寫出該等差數列為：
$$
x, x+1, x+2+,...,x+(k-1)
$$
經過化減後：
\begin{align}
kx&+\frac{(k-1)k}{2}=N \\
kx&=N-\frac{(k-1)k}{2} \\
kx&=N-[1+2+...+(k-1)]
\end{align}

因此我們就可以透過 $\{N-[1+2+...+(k-1)]\}\ mod\ k=0$ 來判斷這個 k 長度的連續整數是否是一組解。

第 5 行先初始化回傳值為 1，因為任何大於 0 的整數都必定會有自己的一個解。
第 7 行的迴圈就是從 k=2 的長度開始找起，第 8 行就是 $\{N-[1+2+...+(k-1)]\}$，所以 `ZZZ` = `x += (1 - k)`，第 9 行就是判斷能不能整除，可以的話就代表可以找到一個解 `ret++`。


### 延伸問題
:::info
TODO
- 嘗試改寫為更有效率的實作
:::

## 參考資料
* [IEEE754 wiki](https://zh.wikipedia.org/wiki/IEEE_754)
* [IEEE754 二進制浮點數種類](https://edisonx.pixnet.net/blog/post/46103946)
* [浮點數運算和定點數操作](https://hackmd.io/@NwaynhhKTK6joWHmoUxc7Q/H1SVhbETQ)
* [以牛頓法求整數開二次方根的近似值 (第 19 頁)](http://science.hsjh.chc.edu.tw/upload_works/106/44f6ce7960aee6ef868ae3896ced46eb.pdf)
* [牛頓求根法](https://www.youtube.com/watch?v=Quw4ZHLH2CY&ab_channel=CUSTCourses)
* [Union type wiki](https://en.wikipedia.org/wiki/Union_type)
* [你所不知道的C語言：技巧篇](https://hackmd.io/@sysprog/c-trick?type=view#do----while0-%E5%B7%A8%E9%9B%86)
* [macro do{}while(0)](http://github.tiankonguse.com/blog/2014/09/30/do-while.html)
* [How to calculate a square root without a calculator](https://www.homeschoolmath.net/teaching/square-root-algorithm.php)
* [Why the square root algorithm works](https://www.homeschoolmath.net/teaching/sqr-algorithm-why-works.php)
* [Square Root of Binary Number Tutorial](https://www.youtube.com/watch?v=G_4HE3ek4c4&ab_channel=ShellWave)
* [Methods of computing square roots](https://en.wikipedia.org/wiki/Methods_of_computing_square_roots#Binary_numeral_system_(base_2))
* [快速整數平方根算法](http://www.azillionmonkeys.com/qed/ulerysqroot.pdf)