mathematics

Python 上的 reals 套件 (需要 3.10+ 以上才能裝)

看到「A lightweight python3 library for arithmetic with real numbers.」這個有趣的 Python 延伸套件，可以用他進行高精度的實數運算...

一開始在 Python 3.9 環境裝，結果就跳出需要 3.10+ 的環境，想了一下，開了一個 Docker container 裝 pyenv 來測，測過以後覺得還蠻有趣的，看起來之後把預設環境變成 3.10+ 應該會裝起來用...

這個 reals 的重點在於保證顯示數字的正確性：

It allows you to compute approximations to an arbitrary degree of precision, and, contrary to most other libraries, guarantees that all digits it displays are correct.

目前支援的常數與操作有這些：

Constants: pi, e, phi
Functions related to powers: sqrt, exp, log
Operators: negation, addition, subtraction, multiplication, division, powers
Trigonometric functions: sin, sinh, csc, csch, cos, cosh, sec, sech, tan, tanh, cot, coth

用法的部份，先把 reals 拉進來：

>>> from reals import sqrt

然後用法算直覺：

>>> sqrt2 = sqrt(2)
>>> sqrt2
<reals._real.Real object at 0x10d182560 (approximate value: 1.41421)>

>>> sqrt2.evaluate(10)
'1.4142135624'
>>> '{:.10f}'.format(sqrt2)
'1.4142135624'
>>> sqrt2.to_decimal(10)
Decimal('1.4142135624')

不過作者有提到效能沒有處理到很好，所以應該是拿來快速做一些運算得到結果而已。

用 GPT-3 解讀程式碼

在 Hacker News 上看到的方法，Simon Willison 試著把程式碼餵進 GPT-3，然後問 GPT-3 程式碼的意思，看起來答的還不錯：「Using GPT-3 to explain how code works」，對應的討論 (包括 Simon Willison 的回應) 則可以在「Using GPT-3 to explain how code works (simonwillison.net)」這邊看到。

第一個範例裡面可以解讀 regular expression，雖然裡面對 (?xm) 的解讀是錯的，但我會說已經很強了...

第二個範例在解釋 Shadow DOM，看起來也解釋的很不錯...

第三個範例回來原來產生程式碼的例子，拿來生 SQL 指令。

後面的 bonus 題目居然是拿來解釋數學公式，他直接丟 TeX 文字進去要 GPT-3 解釋柯西不等式 (Cauchy–Schwarz inequality)。這樣我想到以前高微作業常常會有一堆證明題，好像可以丟進去要 GPT-3 給證明耶...

從三角函數 cosine 的實做問題學一些週邊知識...

前幾天在 Hacker News 上看到「Implementing Cosine in C from Scratch (2020) (austinhenley.com)」這篇 2020 的文章，原文是「Implementing cosine in C from scratch」，裡面內在講自己刻三角函數的 cosine 所遇到的一些嘗試。

cosine 是很基本的函數，所以可以使用的地方很多。另外一方面，也因為他不是那麼直覺就可以實做出來，在現代的實做裡面其實藏了超多細節...

不過真的有趣的是在翻 Hacker News 上的討論時陸陸續續翻其他的資料看到的知識。

第一個看到的是 Intel 對於 FPU-based 指令集內的 FSIN 因為 π 的精度不夠而導致誤差超大 (尤其是在 0 點附近的時候)：「Intel Underestimates Error Bounds by 1.3 quintillion」，然後 AMD 是「相容」到底，所以一樣慘：「Accuracy of FSIN and other x87 trigonometric instructions on AMD processors」。

這個就是有印象，但是太久沒有提到就會忘記...

第二個是 musl libc 裡的 cosine 實做 (看註解應該是從 FreeBSD 的 libc 移植過來的？)：「__cos.c\math\src」與「cos.c\math\src」(話說 cgit 在 html 內 title 的內容對路徑的表達方式頗有趣，居然是反過來放...)。

拆開的部份是先將範圍限制在 $[-\pi/4, \pi/4]$ 後 (這個部份看起來是透過 __rem_pio2.c 處理)，再丟進公式實際運算。

另外帶出來第三個知識，查資料的時候翻到 binary64 (這也是 C 語言裡面的 double) 與 binary128 的差異：

而大家很常拿來惡搞的 double double 則是利用兩個 double 存放，形式是 $v = head + tail$ ，利用不同的 exponent 表示來不同部份的值，以提高經度：

A common software technique to implement nearly quadruple precision using pairs of double-precision values is sometimes called double-double arithmetic.

不過這樣的精確度只能到 106 bits，雖然跟 binary128 能達到的 113 bits 相比低了一些，但在大多數的情況下也還算夠用：

Using pairs of IEEE double-precision values with 53-bit significands, double-double arithmetic provides operations on numbers with significands of at least[4] 2 × 53 = 106 bits (...), only slightly less precise than the 113-bit significand of IEEE binary128 quadruple precision.

Python 裡使用超過 Double Precision 的運算

Hacker News 上爬到的，是一篇 2019 的文章：「When Double Precision Is Not Enough (adambaskerville.github.io)」，原文在「T>T: When Double Precision is Not Enough」。

作者是拿矩陣 (matrix) 的運算當例子，遇到了 double precision 造成的計算誤差問題：

There is no error with the program; this discrepancy is caused by a loss of numerical accuracy in the eigenvalue calculation due to the limitation of hardware double precision (16-digit).

解法是用 mpmath 增加精度，算是一種暴力解，到要注意計算會慢很多：

Note that this library is incredibly slow for large matrices, so is best avoided for most applications.

另外在 Hacker News 的討論串裡面看到個有趣的東西：「Herbie: Automatically Improving Floating Point Accuracy」這個專案，你把公式丟進去，Herbie 會試著提供等價公式來維持精度，像是 $\sqrt{x+1} - \sqrt{x} = 1 / ( \sqrt{x+1} + \sqrt{x} )$ 這種東西。

半自動化幫你改善...

在電腦上計算二進制圓周率 π 的公式

在 Hacker News 的最新列表上面看到的演算法，用來計算圓周率 π 的公式：「Bellard's formula」。

$\pi = \frac1{2^6} \sum_{n=0}^\infty \frac{(-1)^n}{2^{10n}} \, ( -\frac{2^5}{4n+1} - \frac1{4n+3} + \frac{2^8}{10n+1} - \frac{2^6}{10n+3} - \frac{2^2}{10n+5} - \frac{2^2}{10n+7} + \frac1{10n+9} )$

這個公式特別的地方在於 $\frac{1}{2^{10n}}$ 的部份，這使得整個公式可以很迅速的算出某個二進制位置上的值 (要做一些條件判斷，相較於十進制轉二進制的方式快多了)，而且可以馬上想到一些平行運算的方式...

另外一個讓我注意到的點是，這居然是 Fabrice Bellard 丟出來的公式，真的是到處看到他的蹤跡耶...

第 48 個梅森質數確認

在「Great Internet Mersenne Prime Search - PrimeNet」這邊看到 GIMPS 確認了 $M_{57885161}$ (也就是 $2^{57885161}-1$ ) 是第 48 個梅森質數：

All tests smaller than the 48th Mersenne Prime, M(57885161), have been verified

第 47 個梅森質數是 $M_{43112609}$ ，但一直都沒有確認完到 $M_{57885161}$ 中間是否還有梅森質數，這次宣佈全部都算完，而且雙重確認過，所以可以確認第 48 個梅森質數。

下一個要確認的是第 49 個梅森質數，目前已知的下一個是 $M_{74207281}$ ，也就是要看到 $M_{74207281}$ 中間是否還有其他的梅森質數...

各種反直覺的項目

前陣子看到「The most counterintuitive facts in all of mathematics, computer science, and physics」這篇，講各種反直覺的項目，有空的時候拉一兩個來讀還蠻有趣的...

第一個提到的是 Homomorphic encryption。在密碼學的保護概念中，密文是不能被操作的，但我們可以透過重新設計密碼學系統，讓密文可以被運算：

1. It is possible to compute over encrypted data without access to the secret key: https://en.wikipedia.org/wiki/Homomorphic_encryption

2009 年的時候第一個 FHE (fully homomorphic encryption) 演算法被提出來後就被大量討論，這代表你可以把資料加密後丟上雲端，利用雲端的大量運算資源運算完後，再把資料拉回地端解密，得到運算後的結果。這算是 Lattice-based cryptography 的一個很重大的應用。

其他的主題也蠻有趣的，先丟著找機會翻...

在 Hacker News 上看到選擇公理

沒想到會在 Hacker News 的首頁上看到這麼硬核的主題，選擇公理 (Axiom of choice，通常縮寫成 AC)：「What is the Axiom of Choice?」，對應的討論在「What is the Axiom of Choice? (jaydaigle.net)」。

出自「xkcd: Set Theory」

應該是大一教集合論的時候學到的，算是一個非常重要的公設，雖然的確有些數學系統是可以假定 AC 不成立，但用起來會不太好用，主要是因為「對於集合 $S$ ，取出任意一個元素」這類用法太常出現，在沒有 AC 的情況下這件事情就不一定能操作了...

我們目前常用的數學一般是建立在 Zermelo-Fraenkel Set Theory (ZF) 這個公理系統加上 AC，簡寫變成 ZFC。而 AC 在集合論常常會被拿出來說明，主要還是因為在歷史上花了不少力氣才證明 ZF 與 AC 的相對協調性 (ZF 與 AC 不衝突)，以及 ZF 與 AC 獨立性 (ZF 無法推導出 AC)。

有了 AC 後就會再解釋連續統假設 (Continuum hypothesis，簡稱 CH)，也就是 $\mathbb{N}$ 與 $2^{\mathbb{N}}$ 之間存不存在一個集合 $S$ 使得 $|\mathbb{N}| < |S| < |2^{\mathbb{N}}|$ 。