Hamming numbers: Difference between revisions
Line 151: | Line 151: | ||
(cons 1 (smerge3 (map*n 2 hamming) (map*n 3 hamming) (map*n 5 hamming)))))</lang> |
(cons 1 (smerge3 (map*n 2 hamming) (map*n 3 hamming) (map*n 5 hamming)))))</lang> |
||
Note that this version uses a lot of space and time after calculating a few hundred thousand elements of the sequence. This is no doubt due to its "holding on to the head": it maintains the entire generated sequence in memory. |
Note that this version uses a lot of space and time after calculating a few hundred thousand elements of the sequence. This is no doubt due to its "holding on to the head": it maintains the entire generated sequence in memory. |
||
=={{header|D}}== |
|||
D V2 version, from the Java version. This version keeps the whole sequence in memory. |
|||
<lang d>import std.stdio: writeln; |
|||
import std.bigint: BigInt; |
|||
import std.algorithm: min, map; |
|||
import std.range: iota; |
|||
import std.array: array; |
|||
BigInt hamming(int limit) { |
|||
BigInt two = 2; |
|||
BigInt three = 3; |
|||
BigInt five = 5; |
|||
auto h = new BigInt[limit]; |
|||
h[0] = 1; |
|||
BigInt x2 = 2; |
|||
BigInt x3 = 3; |
|||
BigInt x5 = 5; |
|||
int i, j, k; |
|||
foreach (ref el; h[1 .. $]) { |
|||
el = min(x2, x3, x5); |
|||
if (x2 == el) |
|||
x2 = two * h[++i]; |
|||
if (x3 == el) |
|||
x3 = three * h[++j]; |
|||
if (x5 == el) |
|||
x5 = five * h[++k]; |
|||
} |
|||
return h[$ - 1]; |
|||
} |
|||
const(char)[] bigIntRepr(BigInt i) { |
|||
const(char)[] result; |
|||
i.toString((const(char)[] s){ result = s; }, "d"); |
|||
return result; |
|||
} |
|||
void main() { |
|||
writeln(array(map!bigIntRepr(map!hamming(iota(1, 21))))); |
|||
writeln(bigIntRepr(hamming(1691))); |
|||
writeln(bigIntRepr(hamming(1_000_000))); |
|||
}</lang> |
|||
Output: |
|||
<pre>1 2 3 4 5 6 8 9 10 12 15 16 18 20 24 25 27 30 32 36 |
|||
2125764000 |
|||
519312780448388736089589843750000000000000000000000000000000000000000000000000000000</pre> |
|||
=={{header|Factor}}== |
=={{header|Factor}}== |
Revision as of 12:02, 10 June 2010
![Task](http://static.miraheze.org/rosettacodewiki/thumb/b/ba/Rcode-button-task-crushed.png/64px-Rcode-button-task-crushed.png)
You are encouraged to solve this task according to the task description, using any language you may know.
Hamming numbers are numbers of the form
- .
Generate the sequence of Hamming numbers, in increasing order. In particular:
- Show the first twenty Hamming numbers.
- Show the 1691st Hamming number (the last one below ).
- Show the one millionth Hamming number (if the language – or a convenient library – supports arbitrary-precision integers).
References
- wp:Hamming_numbers
- wp:Smooth_number
- Hamming problem from Dr. Dobb's CodeTalk
ALGOL 68
Hamming numbers are generated in a trivial iterative way as in the Python version below. This program keeps the series needed to generate the numbers as short as possible using flexible rows; on the downside, it spends considerable time on garbage collection.
<lang algol68> PR precision=100 PR
MODE SERIES = FLEX [1 : 0] UNT, # Initially, no elements #
UNT = LONG LONG INT; # A 100-digit unsigned integer #
PROC hamming number = (INT n) UNT: # The n-th Hamming number #
CASE n IN 1, 2, 3, 4, 5, 6, 8, 9, 10, 12 # First 10 in a table # OUT # Additional operators # OP MIN = (INT i, j) INT: (i < j | i | j), MIN = (UNT i, j) UNT: (i < j | i | j); PRIO MIN = 9; OP LAST = (SERIES h) UNT: h[UPB h]; # Last element of a series # OP +:= = (REF SERIES s, UNT elem) VOID: # Extend a series by one element, only keep the elements you need # (INT lwb = (i MIN j) MIN k, upb = UPB s; REF SERIES new s = HEAP FLEX [lwb : upb + 1] UNT; (new s[lwb : upb] := s[lwb : upb], new s[upb + 1] := elem); s := new s ); # Determine the n-th hamming number iteratively # SERIES h := 1, # Series, initially one element # UNT m2 := 2, m3 := 3, m5 := 5, # Multipliers # INT i := 1, j := 1, k := 1; # Counters # TO n - 1 DO h +:= (m2 MIN m3) MIN m5; (LAST h = m2 | m2 := 2 * h[i +:= 1]); (LAST h = m3 | m3 := 3 * h[j +:= 1]); (LAST h = m5 | m5 := 5 * h[k +:= 1]) OD; LAST h ESAC;
FOR k TO 20 DO print ((whole (hamming number (k), 0), blank)) OD; print ((newline, whole (hamming number (1 691), 0))); print ((newline, whole (hamming number (1 000 000), 0))) </lang>
Sample output:
1 2 3 4 5 6 8 9 10 12 15 16 18 20 24 25 27 30 32 36 2125764000 519312780448388736089589843750000000000000000000000000000000000000000000000000000000
C
C does not natively support long numbers. The following will print the first 20 and the 1690th hamming numbers <lang c>#include <stdio.h>
typedef unsigned vlong;
typedef void (*Callback)(int,vlong);
void vextend( vlong *list, int *ix, vlong val) {
vlong *val_l = list + *ix-1; vlong *val_h = list + *ix;
while (*val_l > val) --val_l; if (*val_l == val) return;
val_l = list + *ix-1; if (*ix<8000-1) *ix+=1; else if (val >= *val_h) return;
while ( *val_l > val) *val_h-- = *val_l--; *val_h = val;
}
unsigned MAX_UL = 0XFFFFFFFF;
void hamming( int n, Callback cbak) {
vlong explist[8000]; int ix1, ix2; int k2=0,k3=0,k5=0; vlong v; explist[0] = 1; ix1 = 0; ix2 = 1;
while ((ix1 < ix2) && (ix1 < n) ) { v = explist[ix1]; (*cbak)(ix1,v); if ( v < MAX_UL/(unsigned)2) vextend(explist, &ix2, 2*v); if ( v < MAX_UL/(unsigned)3) vextend(explist, &ix2, 3*v); if ( v < MAX_UL/(unsigned)5) vextend(explist, &ix2, 5*v); ix1++; }
}
void hprint(int ix, vlong val) {
if ((ix < 20) || (ix=1689) || (ix==1847)) printf("h[%4d]= %u\n", 1+ix,val);
}
int main(int argc, char *argv[]) {
printf("%u\n", MAX_UL); hamming(5002, &hprint); return 0;
}</lang>
Clojure
This version implements Dijkstra's merge solution, so is closely related to the Haskell version. <lang clojure>(defn smerge [xs ys]
(lazy-seq (let [x (first xs), y (first ys), [z xs* ys*] (cond (< x y) [x (rest xs) ys] (> x y) [y xs (rest ys)] :else [x (rest xs) (rest ys)])] (cons z (smerge xs* ys*)))))
(defn smerge3 [xs ys zs]
(smerge xs (smerge ys zs)))
(defn map*n [n ks] (map #(* n %) ks))
(def hamming
(lazy-seq (cons 1 (smerge3 (map*n 2 hamming) (map*n 3 hamming) (map*n 5 hamming)))))</lang>
Note that this version uses a lot of space and time after calculating a few hundred thousand elements of the sequence. This is no doubt due to its "holding on to the head": it maintains the entire generated sequence in memory.
D
D V2 version, from the Java version. This version keeps the whole sequence in memory. <lang d>import std.stdio: writeln; import std.bigint: BigInt; import std.algorithm: min, map; import std.range: iota; import std.array: array;
BigInt hamming(int limit) {
BigInt two = 2; BigInt three = 3; BigInt five = 5;
auto h = new BigInt[limit]; h[0] = 1; BigInt x2 = 2; BigInt x3 = 3; BigInt x5 = 5; int i, j, k;
foreach (ref el; h[1 .. $]) { el = min(x2, x3, x5); if (x2 == el) x2 = two * h[++i]; if (x3 == el) x3 = three * h[++j]; if (x5 == el) x5 = five * h[++k]; } return h[$ - 1];
}
const(char)[] bigIntRepr(BigInt i) {
const(char)[] result; i.toString((const(char)[] s){ result = s; }, "d"); return result;
}
void main() {
writeln(array(map!bigIntRepr(map!hamming(iota(1, 21))))); writeln(bigIntRepr(hamming(1691))); writeln(bigIntRepr(hamming(1_000_000)));
}</lang> Output:
1 2 3 4 5 6 8 9 10 12 15 16 18 20 24 25 27 30 32 36 2125764000 519312780448388736089589843750000000000000000000000000000000000000000000000000000000
Factor
<lang factor>USING: accessors deques dlists fry kernel make math math.order
IN: rosetta.hamming
TUPLE: hamming-iterator 2s 3s 5s ;
- <hamming-iterator> ( -- hamming-iterator )
hamming-iterator new 1 1dlist >>2s 1 1dlist >>3s 1 1dlist >>5s ;
- enqueue ( n hamming-iterator -- )
[ [ 2 * ] [ 2s>> ] bi* push-back ] [ [ 3 * ] [ 3s>> ] bi* push-back ] [ [ 5 * ] [ 5s>> ] bi* push-back ] 2tri ;
- next ( hamming-iterator -- n )
dup [ 2s>> ] [ 3s>> ] [ 5s>> ] tri 3dup [ peek-front ] tri@ min min [ '[ dup peek-front _ = [ pop-front* ] [ drop ] if ] tri@ ] [ swap enqueue ] [ ] tri ;
- next-n ( hamming-iterator n -- seq )
swap '[ _ [ _ next , ] times ] { } make ;
- nth-from-now ( hamming-iterator n -- m )
1 - over '[ _ next drop ] times next ;</lang>
<hamming-iterator> 20 next-n . <hamming-iterator> 1691 nth-from-now . <hamming-iterator> 1000000 nth-from-now .
Lazy lists aren quite slow in Factor, but still. <lang factor>USING: combinators fry kernel lists lists.lazy locals math ; IN: rosetta.hamming-lazy
- sort-merge ( xs ys -- result )
xs car :> x ys car :> y { { [ x y < ] [ [ x ] [ xs cdr ys sort-merge ] lazy-cons ] } { [ x y > ] [ [ y ] [ ys cdr xs sort-merge ] lazy-cons ] } [ [ x ] [ xs cdr ys cdr sort-merge ] lazy-cons ] } cond ;
- hamming ( -- hamming )
f :> h! [ 1 ] [ h 2 3 5 [ '[ _ * ] lazy-map ] tri-curry@ tri sort-merge sort-merge ] lazy-cons h! h ;</lang>
20 hamming ltake list>array . 1690 hamming lnth . 999999 hamming lnth .
Haskell
We omit the initial 0 for convenience <lang haskell>hamming = 1 : map (2*) hamming `merge` map (3*) hamming `merge` map (5*) hamming
where merge (x:xs) (y:ys) | x < y = x : xs `merge` (y:ys) | x > y = y : (x:xs) `merge` ys | otherwise = x : xs `merge` ys
main = do
print $ take 20 hamming print $ hamming !! 1690 print $ hamming !! 999999</lang>
Icon and Unicon
Icon
This Icon solution uses Unicon extensions. An Icon only version has not been provided.
Unicon
This solution uses lazy evaluation to improve performance.
<lang Unicon>
- Lazily generate the three Hamming numbers that can be derived directly
- from a known Hamming number h
class Triplet : Class (cv, ce)
method nextVal() suspend cv := @ce end
initially (baseNum) cv := 2*baseNum ce := create (3|5)*baseNum
end
- Generate Hamming numbers, in order. Default is first 30
- But an optional argument can be used to generate more (or less)
- e.g. hamming 5000 generates the first 5000.
procedure main(args)
limit := integer(args[1]) | 30 every write("\t", generateHamming() \ limit)
end
- Do the work. Start with known Hamming number 1 and maintain
- a set of triplet Hamming numbers as they get derived from that
- one. Most of the code here is to figure out which Hamming
- number is next in sequence (while removing duplicates)
procedure generateHamming()
triplers := set() insert(triplers, Triplet(1))
suspend 1 repeat { # Pick a Hamming triplet that *may* have the next smallest number t1 := !triplers # any will do to start
every t1 ~=== (t2 := !triplers) do { if t1.cv > t2.cv then { # oops we were wrong, switch assumption t1 := t2 } else if t1.cv = t2.cv then { # t2's value is a duplicate, so # advance triplet t2, if none left in t2, remove it t2.nextVal() | delete(triplers, t2) } }
# Ok, t1 has the next Hamming number, grab it suspend t1.cv insert(triplers, Triplet(t1.cv)) # Advance triplet t1, if none left in t1, remove it t1.nextVal() | delete(triplers, t1) }
end </lang>
J
Solution:
A concise tacit expression using a (right) fold:
<lang j>hamming=: {. (/:~@~.@], 2 3 5 * {)/ @ (1x,~i.@-)</lang>
Example usage: <lang j> hamming 20 1 2 3 4 5 6 8 9 10 12 15 16 18 20 24 25 27 30 32 36
{: hamming 1691
2125764000</lang>
For the millionth (and billionth (1e9)) Hamming number see the hn
verb from Hamming Number essay on the J wiki.
Explanation:
I'll explain this J-sentence by dividing it in three parts from left to right omitting the leftmost {.
:
- sort and remove duplicates
<lang j> /:~@~.@]</lang>
- produce (the next) 3 elements by selection and multiplication:
<lang j> 2 3 5 * {</lang> or <lang j> 2 3 5 * LHA { RHA</lang>
- the RH part forms an array of descending indices and the initial Hamming number 1
<lang j> (1x,~i.@-)</lang> e.g. if we want the first 5 Hamming numbers, it produces the array:
4 3 2 1 0 1
Java
It's recommended to use option -Xmx256m in case of problems calculating the 1,000,000th number.
<lang java> import java.math.BigInteger; public static BigInteger hamming(int limit) {
final BigInteger TWO=BigInteger.valueOf(2); final BigInteger THREE=BigInteger.valueOf(3); final BigInteger FIVE=BigInteger.valueOf(5);
final BigInteger h[]=new BigInteger[limit]; h[0]=BigInteger.ONE; int i=0,j=0,k=0;
BigInteger x2=TWO.multiply(h[i]); BigInteger x3=THREE.multiply(h[j]); BigInteger x5=FIVE.multiply(h[k]);
for(int n=1; n<limit; n++) {
h[n]=x2.min(x3.min(x5));
if(x2.equals(h[n])) x2=TWO.multiply(h[++i]);
if(x3.equals(h[n])) x3=THREE.multiply(h[++j]);
if(x5.equals(h[n])) x5=FIVE.multiply(h[++k]); } return h[limit-1];
} </lang>
JavaScript
This does not calculate the 1,000,000th Hamming number.
Note the use of for (x in obj)
to iterate over the properties of an object, versus for each (y in obj)
to iterate over the values of the properties of an object.
<lang javascript>function hamming() {
var queues = {2: [], 3: [], 5: []}; var base; var next_ham = 1; while (true) { yield next_ham;
for (base in queues) {queues[base].push(next_ham * base)}
next_ham = [ queue[0] for each (queue in queues) ].reduce(function(min, val) { return Math.min(min,val) });
for (base in queues) {if (queues[base][0] == next_ham) queues[base].shift()} }
}
var ham = hamming(); var first20=[], i=1;
for (; i <= 20; i++)
first20.push(ham.next());
print(first20.join(', ')); print('...'); for (; i <= 1690; i++)
ham.next();
print(i + " => " + ham.next());</lang>
output
1, 2, 3, 4, 5, 6, 8, 9, 10, 12, 15, 16, 18, 20, 24, 25, 27, 30, 32, 36 ... 1691 => 2125764000
Logo
<lang logo>to init.ham
; queues make "twos [1] make "threes [1] make "fives [1]
end to next.ham
localmake "ham first :twos if less? first :threes :ham [make "ham first :threes] if less? first :fives :ham [make "ham first :fives]
if equal? :ham first :twos [ignore dequeue "twos] if equal? :ham first :threes [ignore dequeue "threes] if equal? :ham first :fives [ignore dequeue "fives]
queue "twos :ham * 2 queue "threes :ham * 3 queue "fives :ham * 5
output :ham
end
init.ham repeat 20 [print next.ham] repeat 1690-20 [ignore next.ham] print next.ham</lang>
Lua
<lang lua>function hiter()
hammings = {1} prev, vals = {1, 1, 1} index = 1 local function nextv() local n, v = 1, hammings[prev[1]]*2
if hammings[prev[2]]*3 < v then n, v = 2, hammings[prev[2]]*3 end if hammings[prev[3]]*5 < v then n, v = 3, hammings[prev[3]]*5 end prev[n] = prev[n] + 1 if hammings[index] == v then return nextv() end index = index + 1 hammings[index] = v return v
end return nextv
end
j = hiter() for i = 1, 20 do
print(j())
end n, l = 0, 0 while n < 2^31 do n, l = j(), n end print(l)</lang>
MUMPS
<lang MUMPS>Hamming(n) New count,ok,next,number,which For which=2,3,5 Set number=1 For count=1:1:n Do . Set ok=0 Set:count<21 ok=1 Set:count=1691 ok=1 Set:count=n ok=1 . Write:ok !,$Justify(count,5),": ",number . For which=2,3,5 Set next(number*which)=which . Set number=$Order(next("")) . Kill next(number) . Quit Quit Do Hamming(2000)
1: 1 2: 2 3: 3 4: 4 5: 5 6: 6 7: 8 8: 9 9: 10 10: 12 11: 15 12: 16 13: 18 14: 20 15: 24 16: 25 17: 27 18: 30 19: 32 20: 36 1691: 2125764000 2000: 8062156800</lang>
Oz
<lang oz>declare
fun lazy {HammingFun} 1|{FoldL1 [{MultHamming 2} {MultHamming 3} {MultHamming 5}] LMerge} end
Hamming = {HammingFun}
fun {MultHamming N} {LMap Hamming fun {$ X} N*X end} end
fun lazy {LMap Xs F} case Xs of nil then nil [] X|Xr then {F X}|{LMap Xr F} end end
fun lazy {LMerge Xs=X|Xr Ys=Y|Yr} if X < Y then X|{LMerge Xr Ys} elseif X > Y then Y|{LMerge Xs Yr} else X|{LMerge Xr Yr} end end
fun {FoldL1 X|Xr F} {FoldL Xr F X} end
in
{ForAll {List.take Hamming 20} System.showInfo} {System.showInfo {Nth Hamming 1690}} {System.showInfo {Nth Hamming 1000000}}</lang>
PicoLisp
<lang PicoLisp>(de hamming (N)
(let (L (1) H) (do N (for (X L X (cadr X)) # Find smallest result (setq H (car X)) ) (idx 'L H NIL) # Remove it (for I (2 3 5) # Generate next results (idx 'L (* I H) T) ) ) H ) )
(println (make (for N 20 (link (hamming N))))) (println (hamming 1691)) (println (hamming 1000000))</lang> Output:
(1 2 3 4 5 6 8 9 10 12 15 16 18 20 24 25 27 30 32 36) 2125764000 519312780448388736089589843750000000000000000000000000000000000000000000000000000000
Python
<lang python>from itertools import islice
def hamming2():
\ This version is based on a snippet from: http://dobbscodetalk.com/index.php?option=com_content&task=view&id=913&Itemid=85
When expressed in some imaginary pseudo-C with automatic unlimited storage allocation and BIGNUM arithmetics, it can be expressed as: hamming = h where array h; n=0; h[0]=1; i=0; j=0; k=0; x2=2*h[ i ]; x3=3*h[j]; x5=5*h[k]; repeat: h[++n] = min(x2,x3,x5); if (x2==h[n]) { x2=2*h[++i]; } if (x3==h[n]) { x3=3*h[++j]; } if (x5==h[n]) { x5=5*h[++k]; } h = 1 _h=[h] # memoized multipliers = (2, 3, 5) multindeces = [0 for i in multipliers] # index into _h for multipliers multvalues = [x * _h[i] for x,i in zip(multipliers, multindeces)] yield h while True: h = min(multvalues) _h.append(h) for (n,(v,x,i)) in enumerate(zip(multvalues, multipliers, multindeces)): if v == h: i += 1 multindeces[n] = i multvalues[n] = x * _h[i] # cap the memoization mini = min(multindeces) if mini >= 1000: del _h[:mini] multindeces = [i - mini for i in multindeces] # yield h</lang>
Sample output:
>>> list(islice(hamming2(), 20)) [1, 2, 3, 4, 5, 6, 8, 9, 10, 12, 15, 16, 18, 20, 24, 25, 27, 30, 32, 36] >>> list(islice(hamming2(), 1689, 1690)) [2123366400] >>> list(islice(hamming2(), 999999, 1000000)) [519312780448388736089589843750000000000000000000000000000000000000000000000000000000]
Alternate version using "Cyclic Iterators"
The original author is Raymond Hettinger and the code was first published here under the MIT license. This implementation is quite memory efficient. <lang python>from itertools import tee, chain, groupby from heapq import merge
def raymonds_hamming():
# Generate "5-smooth" numbers, also called "Hamming numbers" # or "Regular numbers". See: http://en.wikipedia.org/wiki/Regular_number # Finds solutions to 2**i * 3**j * 5**k for some integers i, j, and k.
def no_repeats(it): for k,g in groupby(it): yield k
def deferred_output(): for i in output: yield i
result, p2, p3, p5 = tee(deferred_output(), 4) m2 = (2*x for x in p2) # multiples of 2 m3 = (3*x for x in p3) # multiples of 3 m5 = (5*x for x in p5) # multiples of 5 merged = merge(m2, m3, m5) combined = chain([1], merged) # prepend a starting point output = no_repeats(combined) # eliminate duplicates
return result</lang>
Results are the same as before.
Alternate version from the Java version
This is the fastest of the Python implementations, it uses a lot of memory. <lang python>import psyco
def hamming(limit):
h = [1] * limit x2, x3, x5 = 2, 3, 5 i = j = k = 0
for n in xrange(1, limit): h[n] = min(x2, x3, x5) if x2 == h[n]: i += 1 x2 = 2 * h[i] if x3 == h[n]: j += 1 x3 = 3 * h[j] if x5 == h[n]: k += 1 x5 = 5 * h[k]
return h[-1]
psyco.bind(hamming) print [hamming(i) for i in xrange(1, 21)] print hamming(1691) print hamming(1000000)</lang>
R
Recursively find the Hamming numbers below . Works for larger limits, however to find the value, one needs guess the correct limit <lang R> hamming=function(hamms,limit) { tmp=hamms for(h in c(2,3,5)) { tmp=c(tmp,h*hamms) } tmp=unique(tmp[tmp<=limit]) if(length(tmp)>length(hamms)) { hamms=hamming(tmp,limit) } hamms } sort(hamming(1,limit=2^31)[-1]) </lang>
Ruby
This method
<lang ruby>require 'generator'
- the Hamming number generator
hamming = Generator.new do |generator|
next_ham = 1 queues = { 2 => [], 3 => [], 5 => [] } loop do generator.yield next_ham [2,3,5].each {|m| queues[m] << (next_ham * m)} next_ham = [2,3,5].collect {|m| queues[m][0]}.min [2,3,5].each {|m| queues[m].shift if queues[m][0] == next_ham} end
end</lang>
This method,
, does not require a library module.
<lang ruby>hamming = Enumerator.new do |yielder|
next_ham = 1 queues = { 2 => [], 3 => [], 5 => [] }
loop do yielder << next_ham # or: yielder.yield(next_ham)
[2,3,5].each {|m| queues[m]<< (next_ham * m)} next_ham = [2,3,5].collect {|m| queues[m][0]}.min [2,3,5].each {|m| queues[m].shift if queues[m][0]== next_ham} end
end</lang>
And the "main" part of the task <lang ruby>start = Time.now
idx = 1 hamming.each do |ham|
case idx when (1..20), 1691 p [idx, ham] when 1_000_000 p [idx, ham] break end idx += 1
end
puts "elapsed: #{Time.now - start} seconds"</lang>
output:
[1, 1] [2, 2] [3, 3] [4, 4] [5, 5] [6, 6] [7, 8] [8, 9] [9, 10] [10, 12] [11, 15] [12, 16] [13, 18] [14, 20] [15, 24] [16, 25] [17, 27] [18, 30] [19, 32] [20, 36] [1691, 2125764000] [1000000, 519312780448388736089589843750000000000000000000000000000000000000000000000000000000] elapsed: 143.96875 seconds
Scala
<lang scala>class Hamming extends Iterator[BigInt] {
import scala.collection.mutable.Queue val qs = Seq.fill(3)(new Queue[BigInt]) def enqueue(n: BigInt) = qs zip Seq(2, 3, 5) foreach { case (q, m) => q enqueue n * m } def next = { val n = qs map (_.head) min; qs foreach { q => if (q.head == n) q.dequeue } enqueue(n) n } def hasNext = true qs foreach (_ enqueue 1)
}</lang>
However, the usage of closures adds a significant amount of time. The code below, though a bit uglier because of the repetitions, is twice as fast: <lang scala>class Hamming extends Iterator[BigInt] {
import scala.collection.mutable.Queue val q2 = new Queue[BigInt] val q3 = new Queue[BigInt] val q5 = new Queue[BigInt] def enqueue(n: BigInt) = { q2 enqueue n * 2 q3 enqueue n * 3 q5 enqueue n * 5 } def next = { val n = q2.head min q3.head min q5.head if (q2.head == n) q2.dequeue if (q3.head == n) q3.dequeue if (q5.head == n) q5.dequeue enqueue(n) n } def hasNext = true List(q2, q3, q5) foreach (_ enqueue 1)
}</lang>
Usage:
scala> new Hamming take 20 toList res87: List[BigInt] = List(1, 2, 3, 4, 5, 6, 8, 9, 10, 12, 15, 16, 18, 20, 24, 25, 27, 30, 32, 36) scala> new Hamming drop 1690 next res88: BigInt = 2125764000 scala> new Hamming drop 999999 next res89: BigInt = 519312780448388736089589843750000000000000000000000000000000000000000000000000000000
Scheme
<lang scheme>(define-syntax lons
(syntax-rules () ((_ lar ldr) (delay (cons lar (delay ldr))))))
(define (lar lons)
(car (force lons)))
(define (ldr lons)
(force (cdr (force lons))))
(define (lap proc . llists)
(lons (apply proc (map lar llists)) (apply lap proc (map ldr llists))))
(define (take n llist)
(if (zero? n) (list) (cons (lar llist) (take (- n 1) (ldr llist)))))
(define (llist-ref n llist)
(if (= n 1) (lar llist) (llist-ref (- n 1) (ldr llist))))
(define (merge llist-1 . llists)
(define (merge-2 llist-1 llist-2) (cond ((null? llist-1) llist-2) ((null? llist-2) llist-1) ((< (lar llist-1) (lar llist-2)) (lons (lar llist-1) (merge-2 (ldr llist-1) llist-2))) ((> (lar llist-1) (lar llist-2)) (lons (lar llist-2) (merge-2 llist-1 (ldr llist-2)))) (else (lons (lar llist-1) (merge-2 (ldr llist-1) (ldr llist-2)))))) (if (null? llists) llist-1 (apply merge (cons (merge-2 llist-1 (car llists)) (cdr llists)))))
(define hamming
(lons 1 (merge (lap (lambda (x) (* x 2)) hamming) (lap (lambda (x) (* x 3)) hamming) (lap (lambda (x) (* x 5)) hamming))))
(display (take 20 hamming)) (newline) (display (llist-ref 1691 hamming)) (newline) (display (llist-ref 1000000 hamming)) (newline)</lang> Output: <lang>(1 2 3 4 5 6 8 9 10 12 15 16 18 20 24 25 27 30 32 36) 2125764000 out of memory</lang>
Smalltalk
This is a straightforward implementation of the pseudocode snippet found in the Python section. Smalltalk supports arbitrary-precision integers, but the implementation is too slow to try it with 1 million.
<lang smalltalk>Object subclass: Hammer [
Hammer class >> hammingNumbers: howMany [ |h i j k x2 x3 x5| h := OrderedCollection new. i := 0. j := 0. k := 0. h add: 1. x2 := 2. x3 := 2. x5 := 5. [ ( h size) < howMany ] whileTrue: [ |m| m := { x2. x3. x5 } sort first. (( h indexOf: m ) = 0) ifTrue: [ h add: m ]. ( x2 = (h last) ) ifTrue: [ i := i + 1. x2 := 2 * (h at: i) ]. ( x3 = (h last) ) ifTrue: [ j := j + 1. x3 := 3 * (h at: j) ]. ( x5 = (h last) ) ifTrue: [ k := k + 1. x5 := 5 * (h at: k) ]. ]. ^ h sort ]
].
(Hammer hammingNumbers: 20) displayNl. (Hammer hammingNumbers: 1690) last displayNl.</lang>
Tcl
This uses coroutines to simplify the description of what's going on.
<lang tcl>package require Tcl 8.6
- Simple helper: Tcl-style list "map"
proc map {varName list script} {
set l {} upvar 1 $varName v foreach v $list {lappend l [uplevel 1 $script]} return $l
}
- The core of a coroutine to compute the product of a hamming sequence.
- Tricky bit: we don't automatically advance to the next value, and instead
- wait to be told that the value has been consumed (i.e., is the result of
- the [yield] operation).
proc ham {key multiplier} {
global hammingCache set i 0 yield [info coroutine] # Cannot use [foreach]; that would take a snapshot of the list in # the hammingCache variable, so missing updates. while 1 {
set n [expr {[lindex $hammingCache($key) $i] * $multiplier}] # If the number selected was ours, we advance to compute the next if {[yield $n] == $n} { incr i }
}
}
- This coroutine computes the hamming sequence given a list of multipliers.
- It uses the [ham] helper from above to generate indivdual multiplied
- sequences. The key into the cache is the list of multipliers.
- Note that it is advisable for the values to be all co-prime wrt each other.
proc hammingCore args {
global hammingCache set hammingCache($args) 1 set hammers [map x $args {coroutine ham$x,$args ham $args $x}] yield while 1 {
set n [lindex $hammingCache($args) [incr i]-1] lappend hammingCache($args) \ [tcl::mathfunc::min {*}[map h $hammers {$h $n}]] yield $n
}
}
- Assemble the pieces so as to compute the classic hamming sequence.
coroutine hamming hammingCore 2 3 5
- Print the first 20 values of the sequence
for {set i 1} {$i <= 20} {incr i} {
puts [format "hamming\[%d\] = %d" $i [hamming]]
} for {} {$i <= 1690} {incr i} {set h [hamming]} puts "hamming{1690} = $h" for {} {$i <= 1000000} {incr i} {set h [hamming]} puts "hamming{1000000} = $h"</lang> Produces this output:
hamming{1} = 1 hamming{2} = 2 hamming{3} = 3 hamming{4} = 4 hamming{5} = 5 hamming{6} = 6 hamming{7} = 8 hamming{8} = 9 hamming{9} = 10 hamming{10} = 12 hamming{11} = 15 hamming{12} = 16 hamming{13} = 18 hamming{14} = 20 hamming{15} = 24 hamming{16} = 25 hamming{17} = 27 hamming{18} = 30 hamming{19} = 32 hamming{20} = 36 hamming{1690} = 2123366400 hamming{1000000} = 519312780448388736089589843750000000000000000000000000000000000000000000000000000000
Ursala
Smooth is defined as a second order function taking a list of primes and returning a function that takes a natural number to the -th smooth number with respect to them. An elegant but inefficient formulation based on the J solution is the following. <lang Ursala>
- import std
- import nat
smooth"p" "n" = ~&z take/"n" nleq-< (rep(length "n") ^Ts/~& product*K0/"p") <1> </lang> This test program <lang Ursala> main = smooth<2,3,5>* nrange(1,20) </lang> yields this list of the first 20 Hamming numbers.
<1,2,3,4,5,6,8,9,10,12,15,16,18,20,24,25,27,30,32,36>
Although all calculations are performed using unlimited precision, the version above is impractical for large numbers. A more hardcore approach is the following. <lang Ursala>#import std
- import nat
smooth"p" "n" =
~&H\"p" *-<1>; @NiXS ~&/(1,1); ~&ll~="n"->lr -+
^\~&rlPrrn2rrm2Zlrrmz3EZYrrm2lNCTrrm2QAX*rhlPNhrnmtPA2XtCD ~&lrPrhl2E?/~&l ^|/successor@l ~&hl, ^|/~& nleq-<&l+ * ^\~&r ~&l|| product@rnmhPX+-
- cast %nL
main = smooth<2,3,5>* nrange(1,20)--<1691,1000000></lang> The program produces this output. The great majority of time is spent calculating the millionth Hamming number.
< 1, 2, 3, 4, 5, 6, 8, 9, 10, 12, 15, 16, 18, 20, 24, 25, 27, 30, 32, 36, 2125764000, 519312780448388736089589843750000000000000000000000000000000000000000000000000000000>