哈明數變體 |
|
哈明數是那些不能被 2、3 和 5 以外的任何質數整除的數。此數列為 Sloane 中的 A051037 (http://www.research.att.com/~njas/sequences/A051037)。前幾個哈明數為
我們的任務是依遞增順序產生這些數。此任務是一個眾所周知的範例,說明了惰性串列為何很有用。事實上,此程式碼和HammingNumbers上的程式碼都使用惰性串列作為解決方案。第三個實作可以在 Mark Jason Dominus 的 Higher Order Perl 中找到 (http://hop.perl.plover.com/)。我要感謝這兩個實作的作者提供此處的構想。
那麼,讓我們看看惰性串列解決方案進行的方式。
我們將惰性串列定義為對串列第一個元素和串列其他元素的一對做出承諾。(大多數其他人定義惰性串列的方式略有不同:將其定義為第一個元素和對串列其他元素承諾的一對。我比較喜歡使用第一種方法,但兩者間沒有太大的差別。)
現在,我們可以古典地像這樣產生哈明數列
hamming := cons(1, merge(map(x -> 2*x, Hamming), merge(map(x -> 3*x, Hamming), map(x -> 5*x, Hamming))));
此處,cons(a, d) 從第一個元素 a 和串列其他元素 d 建立一個串列;map(f, x) 傳回一個串列,其中每個元素都是將函數 f 套用於 x 的各個元素所產生的結果;merge(a, b) 將兩個已排序的串列合併成另一個已排序的串列;x -> n*x 是一個一元函數,會將其輸入乘以 n;而作業會執行充分的惰性化。
此解決方案使用到一個事實,即一個數是哈明數當且僅當它等於 1,或它是哈明數的 2 倍,或它是哈明數的 3 倍,或它是哈明數的 5 倍。
此方法的問題在於此。有些數會透過倍數產生。例如,6 既是哈明數的 2 倍,也是哈明數的 3 倍。因此,要讓此方法發揮作用,您必須以一種方式定義合併,讓其只將共用元素放置到結果串列一次。
然而,還有一個產生哈明數的方法。為此,讓我們將一個數稱為超哈明數,當且僅當其唯一的質數因數為 2 和 3。然後,我們產生 2 的次方數,再從這些次方數產生超哈明數,然後從這些次方數產生哈明數。
2 的次方數很容易:一個數是 2 的次方數當且僅當它等於 1,或它是 2 的次方數的 2 倍。現在請注意,一個數是超哈明數當且僅當它是一個 2 的次方數,或它是超哈明數的 3 倍。類似地,一個數是哈明數當且僅當它是超哈明數,或它是哈明數的 5 倍。請注意,此方法如何只產生每個數列中的每個數一次。
對應的算式如下。
然而這並不能直接執行。這是原因。序列的方程式是正確的,但不足以產生序列,因為你嘗試計算它們時,會得到一個無限迴圈。為什麼會這樣?
第一個序列,two_power 是正確的。但假設你想要取得 very_hamming 列表的第一個元素。為了執行這項工作,你必須計算 two_power 的第一個元素以及 map(x -> 3*x, very_hamming) 的第一個元素,並取兩者較小的元素。但是,要計算 map(x -> 3*x, very_hamming) 的第一個元素,你必須計算 very_hamming 的第一個元素。這是一個無窮遞迴。
我們知道 very_hamming 的第一個元素是 two_power 的第一個元素,1。然而,從定義中無法清楚了解,因為為了得知這一點,你必須知道 map(x -> 3*x, very_hamming) 的第一個元素大於 1。然而,你在知道 very_hamming 的第一個元素前,是無法得知這一點的。
為了解決這項問題,你必須明確地告訴程式 very_hamming 和 hamming 的頭幾個元素。在進行這些變更後,演算法確實可以執行。
現在讓我們看看程式。
第一部分是我從 HammingNumbers 按字面擷取的 bignum 函式庫。第二部分定義 promise 運算,第三部分定義配對。這兩個分別被視為抽象資料類型。第四部分定義一些惰性列表函式,其中包括 map 和 merge。最後一個區塊定義 Hamming 數序列並列印其中一些數。
-- hamming.lua hamming numbers example -- see https://lua-users.dev.org.tw/wiki/HammingNumbers -- this code is unoptimized -- bignums ----------------------------------------------------------- -- bignum library do -- very limited bignum stuff; just enough for the examples here. -- Please feel free to improve it. local base = 1e15 local fmt = "%015.0f" local meta = {} function meta:__lt(other) if self.n ~= other.n then return self.n < other.n end for i = 1, self.n do if self[i] ~= other[i] then return self[i] < other[i] end end end function meta:__eq(other) if self.n == other.n then for i = 1, self.n do if self[i] ~= other[i] then return false end end return true end end function meta:__mul(k) -- If the base where a multiple of all possible multipliers, then -- we could figure out the length of the result directly from the -- first "digit". On the other hand, printing the numbers would be -- difficult. So we accept the occasional overflow. local offset = 0 if self[1] * k >= base then offset = 1 end local carry = 0 local result = {} local n = self.n for i = n, 1, -1 do local tmp = self[i] * k + carry local digit = tmp % base carry = (tmp - digit) / base result[offset + i] = digit end if carry ~= 0 then n = n + 1 if offset == 0 then for i = n, 2, -1 do result[i] = result[i - 1] end end result[1] = carry end result.n = n return setmetatable(result, meta) end -- Added so that Fib would work; other must be a Bignum function meta:__add(other) local result = {} if self.n < other.n then self, other = other, self end local n, m = self.n, other.n local diff = n - m local carry = 0 local result = {} for i = m, 1, -1 do local tmp = self[i + diff] + other[i] + carry if tmp < base then carry = 0 else tmp = tmp - base carry = 1 end result[i + diff] = tmp end for i = diff, 1, -1 do local tmp = self[i] + carry if tmp < base then carry = 0 else tmp = tmp - base carry = 1 end result[i] = tmp end if carry > 0 then n = n + 1 for i = n, 2, -1 do result[i] = result[i - 1] end result[1] = carry end result.n = n return setmetatable(result, meta) end function meta:__tostring() local tmp = {} tmp[1] = ("%.0f"):format(self[1]) for i = 2, self.n do tmp[i] = fmt:format(self[i]) end return table.concat(tmp) end function Bignum(k) return setmetatable({k, n = 1}, meta) end end -- promises ---------------------------------------------------------- function delay(f) return {promise_tag = true; promise_func = f}; end; function force(p) if not p.promise_tag then error("forcing a non-promise object of type "..type(p)) end; local f = p.promise_func; if f then local x = f(); p.promise_val, p.promise_func = x, nil; return x; else return p.promise_val; end; end; -- pairs ------------------------------------------------------------- -- we need only infinite lists so we won't mess with nulls here function cons(a, d) return {pair_car = a, pair_cdr = d}; end; function car(c) return c.pair_car; end; function cdr(c) return c.pair_cdr; end; -- lazy lists -------------------------------------------------------- -- a lazy list is a promise to a pair whose cdr is a lazy list -- access fields (thus forcing the list) function lcar(l) return car(force(l)); end; function lcdr(l) return cdr(force(l)); end; -- map a lazy list function lmap(f, l) return delay(function() return cons(f(lcar(l)), lmap(f, lcdr(l))); end); end; -- merge two nondecreasing lazy lists to a lazy list function lmerge(a, b) return delay(function() local x, y = lcar(a), lcar(b); if x <= y then return cons(x, lmerge(lcdr(a), b)); else return cons(y, lmerge(a, lcdr(b))); end; end); end; -- iterate a lazy list function lforeach(f, l) f(lcar(l)); return lforeach(f, lcdr(l)); end; function lany(f, l) x = f(lcar(l)); if x then return x; else return lany(f, lcdr(l)); end; end; -- main ------------------------------------------------------------ --dofile "dump.lua"; -- sequence of natural numbers local N; N = delay(function() return cons(1, lmap(function(x) return x + 1; end, N)) end); -- sequence of Hamming numbers local timeser = function(x) return function(y) return y * x; end; end; local H2; H2 = delay(function() return cons(Bignum(1), lmap(timeser(2), H2)); end); local H3; H3 = delay(function() return cons(Bignum(1), delay(function() return cons(Bignum(2), lcdr(lcdr(lmerge(lmap(timeser(3), H3), H2))) ); end) ); end); local H5; H5 = delay(function() return cons(Bignum(1), delay(function() return cons(Bignum(2), lcdr(lcdr(lmerge(lmap(timeser(5), H5), H3))) ); end) ); end); local n, m = 1, 500005; lany(function(a) if 0 == n % 50000 or n <= 200 then print(n, a); end; n = n + 1; return m <= n; end, H5); -- END