## Q: Assume x[i] is a single floating pointing number (“float”) complying with IEEE 754. Calculate the geomean of x[i] without FPU.
Ans.
int geomean(float *x, int n) // n is the amount of x[i] {}
$$
\text{GM} = \sqrt[n]{x_1 \cdot x_2 \cdots x_n}
$$
### 解題想法
我的想法是,將輸入的浮點數轉換成定點數,使用定點數進行所有乘法運算,計算整體乘積。最後使用 牛頓法 求出乘積的 $n$ 次方根,得到近似的幾何平均值。
在計算 geomean 之前,需先檢查乘積 $P=\prod x_i$。若 $P<0$ 且 $n$ 為偶數,則實數域中不存在對應的 $n$ 次方根,因此無實數解;只有在 $n$ 為奇數時才可得到實數結果。
```c
int is_geomean_defined(float *x, int n)
{
int negative_count = 0;
int has_zero = 0;
for (int i = 0; i < n; i++) {
uint32_t bits;
memcpy(&bits, &x[i], sizeof(float));
if ((bits >> 31) == 1) {
negative_count++;
} else if (x[i] == 0.0f) {
has_zero = 1;
break;
}
}
if (has_zero) return 2;
if (negative_count & 1)
return (n & 1);
return 1;
}
```
以下是我自己設計的 test case :
* 一般情形:例如 Case 1,數值皆為正,計算幾何平均的標準情境。
* 負數乘積為正:如 Case 2,負數個數為偶數,結果仍為正數,幾何平均合法。
* 負數乘積為負,且開偶數次方根:如 Case 3,結果為虛數,在實數下不合法。
* 負數乘積為負,但開奇數次方根:如 Case 4,結果為負實數,合法。
* 乘積為 0:如 Case 5,包含 0,定義上 GM 為 0,視為合法,且掠過計算幾何平均。
* 整數完美開根測試:如 Case 6,用來檢驗演算法是否能精確算出 GM。
* 一般浮點數誤差測試:如 Case 7、8,對應實際運算中的誤差情境。
```c
float x[] = {1.0f, 2.0f, 3.0f, 4.0f};
float x2[] = {-2.0f, -3.0f};
float x3[] = {-2.0f, 3.0f};
float x4[] = {-2.0f, 3.0f, 1.0f};
float x5[] = {0.0f, 3.0f, 2.0f};
float x6[] = {2.0f, 8.0f};
float x7[] = {1.5f, 2.0f, 3.0f, 4.0f, -1.57f};
float x8[] = {0.5f, 0.5f, 0.5f};
```
先計算有 FPU 的情況下的正確答案:
```
Case 1: Geomean ≈ 2.213364 (FPU)
Case 2: Geomean ≈ 2.449490 (FPU)
Case 3: Geomean undefined (invalid negatives)
Case 4: Geomean ≈ -1.817121 (FPU)
Case 5: Geomean = 0.0 (contains zero)
Case 6: Geomean ≈ 4.000000 (FPU)
Case 7: Geomean ≈ -2.240993 (FPU)
Case 8: Geomean ≈ 0.500000 (FPU)
```
Case 3 是負數且開偶數次方根,所以不能做 geomean.
### 整數版本的幾何平均計算(無 FPU)
在不具浮點硬體支援的環境下,無法使用 float 或 double 型別進行運算。因此若直接將浮點數近似為整數後處理,會犧牲小數精度,導致結果與實際值偏差。
以下為整數版的幾何平均計算流程,使用位元操作解析 IEEE 754 格式,取出 exponent 部分來計算每個值。乘積計算後,再透過牛頓法開 $n$ 次方根。
Reference: [牛頓法](https://hackmd.io/@TonyLinX/r1QzzVbmgl)
```c
int int_X_pow(int X, int n) {
int result = 1;
for (int i = 0; i < n; i++) {
result *= X;
}
if (result==0) {
return 1; // 避免除以0
}
return result;
}
int geomean(float *x, int n)
{
int product = 1;
int negative_count = 0;
int eps = 1;
int MAX_ITR = 30;
for (int i = 0; i < n; i++){
uint32_t bits;
int xi = 1;
memcpy(&bits, &x[i], sizeof(float));
int exponent = ((bits >> 23) & 0xFF) - 127;
int sign = (bits >> 31) & 1;
if (sign) negative_count++;
if (exponent < 0)
xi = 1;
else
xi = 1 << exponent;
product *= xi;
}
int X = product; // 初始值
for (int i = 0; i < MAX_ITR; i++) {
int X_next = ((n - 1) * X + (product/ int_X_pow(X, n-1) )) / n;
if(llabs(X - X_next) < eps) {
break; // 如果 X 已經收斂,則退出迴圈
}
X = X_next;
}
if (negative_count & 1) X = -X;
return X;
}
```
實驗結果(與 FPU 計算對比) :
```
Case 1: Geomean ≈ 2.213364 (FPU), 2 (Int), error: 9.639802 %
Case 2: Geomean ≈ 2.449490 (FPU), 2 (Int), error: 18.350346 %
Case 3: Geomean undefined (invalid negatives)
Case 4: Geomean ≈ -1.817121 (FPU), -1 (Int), error: 44.967880 %
Case 5: Geomean = 0.0 (contains zero)
Case 6: Geomean ≈ 4.000000 (FPU), 4 (Int), error: 0.000000 %
Case 7: Geomean ≈ -2.240993 (FPU), -3 (Int), error: 33.869240 %
Case 8: Geomean ≈ 0.500000 (FPU), 1 (Int), error: 100.000000 %
```
可以觀察到,若僅使用整數運算處理幾何平均,一旦實際結果無法開出整數根,都會有誤差。
### 定點數版本的幾何平均計算(無 FPU)
#### 定點數版本需調整以下兩個部分:
在牛頓法中,每次迭代會用到 $X^{k-1}_n$。與整數版本不同,我們需要確保每次乘法後仍保持在 Q8.8 格式。這可以透過在每次乘積後進行右移操作來完成:
```diff
-int int_X_pow(int X, int n) {
- int result = 1;
+int64_t int_X_pow_fixed_q8_8(int64_t X, int n) {
+
+ int64_t result = 1LL << 8;
+
for (int i = 0; i < n; i++) {
- result *= X;
+ result = (result * X) >> 8;
}
- if (result==0) {
- return 1; // 避免除以0
- }
-
+
return result;
-}
```
整數版本僅透過 1 << exponent 取得數值的整數,未考慮 $1.\text{mantissa}$ 的小數部分,因此在浮點數轉換時會產生明顯誤差。
相對地,在定點數版本中,將 $1.\text{mantissa}$ 精確地轉換為 Q8.8 格式,可以保留了小數的資訊。這樣的處理方式可顯著降低誤差,使整體結果更接近真實浮點數計算。
```diff
-int geomean(float *x, int n)
+float geomean_fixed_q8_8(float *x, int n)
{
- int product = 1;
+ int64_t product = 1LL << 8;
int negative_count = 0;
- int eps = 1;
- int MAX_ITR = 30;
for (int i = 0; i < n; i++){
uint32_t bits;
- int xi = 1;
memcpy(&bits, &x[i], sizeof(float));
- int exponent = ((bits >> 23) & 0xFF) - 127;
- int sign = (bits >> 31) & 1;
+ uint8_t sign = (bits >> 31) & 1;
+ uint8_t exponent = (bits >>23) & 0xFF;
+ uint32_t mantissa = bits & 0x7FFFFF;
+ uint32_t value = (1<<8) + (mantissa>>15);
- if (sign) negative_count++;
- if (exponent < 0)
- xi = 1;
- else
- xi = 1 << exponent;
-
- product *= xi;
- }
- int X = product; // 初始值
- for (int i = 0; i < MAX_ITR; i++) {
- int X_next = ((n - 1) * X + (product/ int_X_pow(X, n-1) )) / n;
- if(llabs(X - X_next) < eps) {
- break;
+ if (sign) negative_count++;
+ if (expnent > 127) {
+ value <<= (expnent - 127);
+ } else if (expnent < 127) {
+ value >>= (127 - expnent);
}
- X = X_next;
+ product = (product * value) >> 8;
}
- if (negative_count & 1) X = -X;
- return X;
+ int64_t X = product; // 初始值
+ const uint8_t MAX_ITR = 30;
+ const int32_t EPS = 1;
+
+ for(int i = 0; i < MAX_ITR; i++){
+ int64_t X_next = ((n - 1) * X + ((product<<8) / int_X_pow_fixed_q8_8(X, n-1))) / n;
+ if(llabs(X_next - X) < EPS) {
+ X = X_next;
+ break;
+ }
+ X = X_next;
+ }
+
+ if (negative_count & 1)
+ return -(X / 256.0f);
+ else
+ return X/256.0f;
}
```
實驗結果(與 FPU 計算對比) :
```
Case 1: Geomean ≈ 2.213364 (FPU), 2 (Int), error: 9.639802 %, Fixed: 2.210938, error_fixed: 0.109624 %
Case 2: Geomean ≈ 2.449490 (FPU), 2 (Int), error: 18.350346 %, Fixed: 2.449219, error_fixed: 0.011067 %
Case 3: Geomean undefined (invalid negatives)
Case 4: Geomean ≈ -1.817121 (FPU), -1 (Int), error: 44.967880 %, Fixed: -1.816406, error_fixed: 0.039316 %
Case 5: Geomean = 0.0 (contains zero)
Case 6: Geomean ≈ 4.000000 (FPU), 4 (Int), error: 0.000000 %, Fixed: 4.000000, error_fixed: 0.000000 %
Case 7: Geomean ≈ -2.240993 (FPU), -3 (Int), error: 33.869240 %, Fixed: -2.238281, error_fixed: 0.120997 %
Case 8: Geomean ≈ 0.500000 (FPU), 1 (Int), error: 100.000000 %, Fixed: 0.500000, error_fixed: 0.000000 %
```
可以觀察到,在使用 Q8.8 格式的定點數實作後,除了整數根的 Case6 之外,其他 Case 都有顯著的改善,誤差都低於 1%,有些甚至達到 0.01% 以下。相較於整數版本因無法保留小數部分導致大幅誤差,定點數能保有部分小數精度,使乘積與開根操作更加穩定且接近真實值。

### 討論 overflow
之前定點數的版本,我是先將所有浮點數轉成 Q8.8 的定點數後,全部相乘在一起得到乘積 $prdouct$ 再去做 geomean 運算,但沒有在每一次相乘過程中限制 $product$ 的大小,所以當 $n$ 很大或 $x[i]$ 很大的時候,在計算 $product$ 的過程就很容易 overflow。
這是我新增的 test case 用來測試方法是否會 overflow :
* $x_{9}$ 與 $x_{10}$ 是模擬 array 裡的浮點數很大的情況。
* $x_{11}$ 與 $x_{12}$ 是模擬 n 很大的情況。
* $x_{13}$, $x_{14}$, $x_{15}$ 是模擬 $x[i]$ 不都是 $2$ 的冪次以及 n 越來越大的情況下,是否能不 overflow 且保持小數點第六位的精度。
```c
static float x9[] = {1e8f, 1e8f, 1e8f};
static float x10[] = {512.0f, 512.0f, 512.0f, 512.0f, 512.0f, 512.0f, 512.0f};
static float x11[100];
static float x12[1000];
static float x13[10];
static float x14[100];
static float x15[512];
void init_overflow_arrays (void) {
for (int i = 0; i < 100; i++) x11[i] = 32.0f;
for (int i = 0; i < 1000; i++) x12[i] = 2.0f;
for (int i = 0; i < 10; i++) x13[i] = 1.0f + i;
for (int i = 0; i < 100; i++) x14[i] = 256.0f + i;
for (int i = 0; i < 512; i++) x15[i] = 256.0f + i;
}
```
#### 實際遇到 overflow 的部份
以下想要討論,哪裡會遇到 overflow 以及具體我改善的想法。
##### 計算 product 的部份
$$
\text{GM} = \sqrt[n]{x_1 \cdot x_2 \cdots x_n}
= \sqrt[n]{product}
= \sqrt[n]{P \times 2^{E}}
= P^{\frac{1}{n}} \times 2^{\frac{E}{n}}
= P^{\frac{1}{n}} \times (2^{a} \times 2^{\frac{r}{n}})
$$
不再等到計算完所有乘積的結果後,才做幾何平均,而是將原本的 $product$ 分成了 $P$ 跟 $2^{E}$ 分開進行幾何運算。具體來說,在計算 $proudct$ 的過程會多兩次檢查,第一次檢查是檢查 $x[i]$ 是否落在我們控制的範圍內,如果超過,就右移 1 bit 以及 `E++`,這樣可以確保不會因為中途有超大的 $x[i]$,導致 overflow。而第二次的檢查是每一次 $P$ 乘完 $x[i]$ 之後,會檢查 $P$ 是否落在我們控制的範圍內,如果超過,就右移 1 bit 以及 `E++`,這樣可以確保 $P$ 不會越來越大,所以最後會得到一個值為 $P \times 2^{E}$,再拿這個值去做牛頓開 n 次方根 。
當你的 $n$ 或 $x[i]$ 夠大的情況下,你的 $E$ 也會很大,所以我們無法直接去做 `1 << E` ,於是拆分成整數的部份 $2^{a}$ 以及小數的部份 $2^{\frac{r}{n}}$,例如:$2^{\frac{513}{512}} = 2^{1} + 2^{\frac{1}{512}}$。
```
Case 9: Value: 4093640704, Product: 4093640704
Value: 4093640704, Product: -6597069766656000
Value: 4093640704, Product: 0
Case 10: Value[0]: 131072, Product: 131072
Value[1]: 131072, Product: 67108864
Value[2]: 131072, Product: 34359738368
Value[3]: 131072, Product: 17592186044416
Value[4]: 131072, Product: 9007199254740992
Value[5]: 131072, Product: 0
Value[6]: 131072, Product: 0
Case 11: Value[0]: 8192, Product: 8192
Value[1]: 8192, Product: 262144
Value[2]: 8192, Product: 8388608
Value[3]: 8192, Product: 268435456
Value[4]: 8192, Product: 8589934592
Value[5]: 8192, Product: 274877906944
Value[6]: 8192, Product: 8796093022208
Value[7]: 8192, Product: 281474976710656
Value[8]: 8192, Product: 9007199254740992
Value[9]: 8192, Product: 0
```
##### 計算 $2^{\frac{r}{n}}$ 的部份
理論上,可以直接透過 `1 << r` 得到根號內的數值,再去做牛頓法開 $n$ 次方根,但在 Case 15 ,會有 overflow 的問題。例如:在 Case 15 中,$n$ 會是 512 ,而 $r$ 是 463,`1 << 463` 會超過 int64_t 可表達的範圍。可以先透過直式除法將 ${\frac{r}{n}}$ 拆成許多可表達的數值,例如:$2^{\frac{7}{8}} = 2^{\frac{1}{2}} \times 2^{\frac{1}{4}} \times 2^{\frac{1}{8}}$,轉成二進制會變成 $0.111_{2}$,將這些許多可表達的數字相乘可以得到近似於 $2^{\frac{r}{n}}$ 的值。
$$
\frac{r}{n} = 0.b_1b_2b_3\dots_2=\sum \frac{b_i}{2^i} , b_i\in\{0,1\}
$$
先看過去的結果,可以觀察到過去的寫法都會 overflow,無法完成 Case9 ~ Case12 :
```
Case 1: Geomean ≈ 2.213364 (FPU), 2 (Int), error: 9.639802 %, Fixed: 2.210938, error_fixed: 0.109624 %
Case 2: Geomean ≈ 2.449490 (FPU), 2 (Int), error: 18.350346 %, Fixed: 2.449219, error_fixed: 0.011067 %
Case 3: Geomean undefined (invalid negatives)
Case 4: Geomean ≈ -1.817121 (FPU), -1 (Int), error: 44.967876 %, Fixed: -1.816406, error_fixed: 0.039310 %
Case 5: Geomean = 0.0 (contains zero)
Case 6: Geomean ≈ 4.000000 (FPU), 4 (Int), error: 0.000000 %, Fixed: 4.000000, error_fixed: 0.000000 %
Case 7: Geomean ≈ -2.240993 (FPU), -3 (Int), error: 33.869255 %, Fixed: -2.238281, error_fixed: 0.120986 %
Case 8: Geomean ≈ 0.500000 (FPU), 1 (Int), error: 100.000000 %, Fixed: 0.500000, error_fixed: 0.000000 %
Case 9: Geomean ≈ 100000000.000000 (FPU), 0 (Int), error: 100.000000 %, Fixed: 0.000000, error_fixed: 100.000000 %
Case 10: Geomean ≈ 512.000000 (FPU), 0 (Int), error: 100.000000 %, Fixed: 0.000000, error_fixed: 100.000000 %
Case 11: Geomean ≈ 32.000000 (FPU), 0 (Int), error: 100.000000 %, Fixed: 0.691406, error_fixed: 97.839355 %
Case 12: Geomean ≈ 2.000000 (FPU), 0 (Int), error: 100.000000 %, Fixed: 0.882812, error_fixed: 55.859375 %
Case 13: Geomean ≈ 4.528728 (FPU), 524288 (Int), error: 11576836.000000 %, Fixed: 17092.140625, error_fixed: 377315.875000 %
Case 14: Geomean ≈ 304.128235 (FPU), 0 (Int), error: 100.000000 %, Fixed: 0.691406, error_fixed: 99.772659 %
Case 15: Geomean ≈ 488.833832 (FPU), 0 (Int), error: 100.000000 %, Fixed: 0.882812, error_fixed: 99.819405 %
```
#### 具體修改的地方
這裡是討論實際實做的內容
##### 修改 product 的部份
$$
\text{GM} = \sqrt[n]{x_1 \cdot x_2 \cdots x_n}
= \sqrt[n]{product}
= \sqrt[n]{P \times 2^{E}}
= P^{\frac{1}{n}} \times 2^{\frac{E}{n}}
= P^{\frac{1}{n}} \times (2^{a} \times 2^{\frac{r}{n}})
$$
過去定點數的版本,會先計算出最後的乘積 $product$,再去做牛頓法開 n 次方根,但這樣容易 overflow,所以作了以下調整 :
* 將 $prodcut$ 分成 $P \times 2^{E}$
* 檢查 $x[i]$ `(value)` 是否落在我們控制的範圍內 `(Q8_SCALE_FACTOR * 2)`,避免突然的超大值,造成 overflow。
* 相乘完得到新的 $P$,會再檢查新的 $P$ 是否落在控制範圍內 `(63 - Q8_SCALE_FACTOR_BITS - 1)`,確保乘上 $x[i+1]$ 不會 overflow。
```diff
+#define SCALE_BITS 8
+#define SCALE_FACTOR (1 << SCALE_BITS)
-float geomean_fixed_q8_8(float *x, int n)
+float geomean_precise(float *x, int n)
{
- int64_t product = 1LL << 8;
+ int64_t P = 1LL << SCALE_BITS; // P 代表是 product
+ int64_t R= 0;
+
+ int E = 0;
int negative_count = 0;
+
for (int i = 0; i < n; i++){
uint32_t bits;
uint8_t sign = (bits >> 31) & 1;
uint8_t expnent = (bits >>23) & 0xFF;
uint32_t mantissa = bits & 0x7FFFFF;
- uint32_t value = (1<<8) + (mantissa>>15);
-
+ uint32_t value = (1 << SCALE_BITS) + (mantissa >> ((SCALE_BITS < 23) ? (23 - SCALE_BITS) : 0));
+ int32_t val = value;
+ int64_t Pb = P;
+
+
if (sign) negative_count++;
- if (expnent > 127) {
- value <<= (expnent - 127);
- } else if (expnent < 127) {
- value >>= (127 - expnent);
+ E += expnent - 127;
+
+ while(value >= (SCALE_FACTOR * 2)) {
+ value >>= 1;
+ E++;
}
- product = (product * value) >> 8;
-
- }
+ P = (P * value) >> SCALE_BITS;
+
+ while(P >= (1LL << (63 - SCALE_BITS - 2))) {
+ P >>= 1;
+ E += 1;
+ }
+
+ }
```
##### 修改 $2^{\frac{r}{n}}$ 的部份
計算完所有的 $x[i]$ 後,將 $P$ 與 $2^{E}$ 分開做牛頓法得到 $n$ 次方根。在計算 $2^{E}$ 的牛頓法,由於在 $n$ 很大的情況下,$E$ 也很大,不能直接 `1 << E`,所以將 $E$ 的部分拆成整數的部份 $2^{a}$ 以及小數的部分 $2^{\frac{r}{n}}$,例如 : $2^{\frac{975}{512}} = 2^{1} \times 2^{\frac{463}{512}}$。整數的部分只需 shift bit 即可完成,但小數的部份須透過直式除法將 ${\frac{r}{n}}$ 拆成許多可表達的數值,例如:$2^{\frac{7}{8}} = 2^{\frac{1}{2}} \times 2^{\frac{1}{4}} \times 2^{\frac{1}{8}}$,轉成二進制會變成 $0.111_{2}$,將這些許多可表達的數字相乘可以得到近似於 $2^{\frac{r}{n}}$ 的值。那這些數值 $2^{\frac{1}{2}}$, $2^{\frac{1}{4}}$, $2^{\frac{1}{8}}$ 可以先建立表格存起來,在計算過程只要透過查表即可完成。
```diff
+int64_t X = 1LL;
+int64_t P_X; //product 的初始值 x
+int64_t E_X = 1; // exponent 的初始值 x
+int64_t a = E / n;
+int64_t r = E % n;
+int64_t X_r = 0;
+int count = r/20;
-int64_t X;
+P_X = newton_precise(P, n);
-X = newton_q8_8(product, n);
+if (r) {
+ E_X = pow2_frac_q23(r, n, 18, SCALE_BITS);
+ X = (P_X * E_X);
+ X >>= SCALE_BITS;
+}
+else
+ X = P_X;
+
+if (a > 0) X = X << a;
+else X = X >> -a;
if (negative_count & 1) X = -X;
-return X / 256.0f;
+return (X) / (float)SCALE_FACTOR;
```
以下是我建立的表格 :
```c
/* ---- Q24.8 (scale = 1<<8) ---- */
static const int64_t k_tbl_8[18] = {
362, 304, 279, 267,
262, 259, 257, 257,
256, 256, 256, 256,
256, 256, 256, 256,
256, 256
};
int64_t pow2_frac_q23(uint32_t p, uint32_t q, int B, int SCALE_BIT)
{
int64_t res = 1 << SCALE_BIT;
for (int i = 0; i < B; ++i) {
p <<= 1;
if (p >= q) {
p -= q;
if (SCALE_BIT == 23)
res = (res * k_tbl_23[i]) >> SCALE_BIT;
else if (SCALE_BIT == 16)
res = (res * k_tbl_16[i]) >> SCALE_BIT;
else if (SCALE_BIT == 8)
res = (res * k_tbl_16[i]) >> SCALE_BIT;
else if (SCALE_BIT == 4)
res = (res * k_tbl_16[i]) >> SCALE_BIT;
}
if (p == 0) break;
}
return res;
}
```
實驗結果 :
```
Case 1: Geomean ≈ 2.213364 (FPU), 2 (Int), Precise: 2.2109375, error_precise: 0.1096243 %
Case 2: Geomean ≈ 2.449490 (FPU), 2 (Int), Precise: -1.8046875, error_precise: 0.6842172 %
Case 5: Geomean = 0.0 (contains zero)
Case 4: Geomean ≈ -1.817121 (FPU), -1 (Int), Precise: -1.8171167, error_precise: 0.0002099 %
Case 6: Geomean ≈ 4.000000 (FPU), 4 (Int), Precise: 4.0000000, error_precise: 0.0000000 %
Case 7: Geomean ≈ -2.240993 (FPU), -3 (Int), Precise: -2.2265625, error_precise: 0.6439131 %
Case 8: Geomean ≈ 0.500000 (FPU), 1 (Int), Precise: 0.5000000, error_precise: 0.0000000 %
Case 9: Geomean ≈ 100000000.000000 (FPU), 0 (Int), Precise: 99876864.0000000, error_precise: 0.1231360 %
Case 10: Geomean ≈ 512.000000 (FPU), 0 (Int), Precise: 512.0000000, error_precise: 0.0000000 %
Case 11: Geomean ≈ 32.000000 (FPU), 0 (Int), Precise: 32.0000000, error_precise: 0.0000000 %
Case 12: Geomean ≈ 2.000000 (FPU), 0 (Int), Precise: 2.0000000, error_precise: 0.0000000 %
Case 13: Geomean ≈ 4.528728 (FPU), 524288 (Int), Precise: 4.4843750, error_precise: 0.9793805 %
Case 14: Geomean ≈ 304.128235 (FPU), 0 (Int), Precise: 303.0000000, error_precise: 0.3709734 %
Case 15: Geomean ≈ 488.833832 (FPU), 0 (Int), Precise: 487.0000000, error_precise: 0.3751442 %
```
由此圖可以觀察到,定點數的版本 (Fixed) 在 Case 1 ~ Case 8 表現很好,但在 Case 9 ~ Case 15 都會 oveflow。經過上述的改善後 (Precise),不但不會 overflow,而且誤差律都可以降到 $1\%$ 以下。

關於精細度的問題可以使用更多位數的定點數使結果更佳精準,例如 Q16 或 Q23。可以觀察到使用 Q23 的誤差率遠低於 Q8。
