Hamming numbers: Difference between revisions

Content deleted Content added

Inline

Revision as of 15:02, 20 December 2010

Hamming numbers are numbers of the form

H=2^{i}\cdot 3^{j}\cdot 5^{k},\;\mathrm {where} \;i,j,k\geq 0

.

Hamming numbers are also known as ugly numbers and also 5-smooth numbers (numbers whose prime divisors are less or equal to 5).

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 $2^{31}$ ).
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

Ada

Works with: GNAT

GNAT provides the datatypes Integer, Long_Integer and Long_Long_Integer.

Values for GNAT Pro 6.3.1, 64 bit Linux version:

Integer covers the range -2**31 .. 2**31-1 (-2147483648 .. 2147483647).
Long_Integer and Long_Long_Integer cover the range -2**63 .. 2**63-1 (-9223372036854775808 .. 9223372036854775807).

Using your own modular integer type (for example type My_Unsigned_Integer is mod 2**64;), you can expand the range to 0 .. 18446744073709551615, but this still is not enough for the millionth Hamming number.

For bigger numbers, you have to use an external library, for example GNU GMP.

The code for calculating the Hamming numbers is kept generic, to easily expand the range by changing the concrete type.

<lang Ada>with Ada.Text_IO; procedure Hamming is

  generic
     type Int_Type is private;
     Zero  : Int_Type;
     One   : Int_Type;
     Two   : Int_Type;
     Three : Int_Type;
     Five  : Int_Type;
     with function "mod" (Left, Right : Int_Type) return Int_Type is <>;
     with function "/"   (Left, Right : Int_Type) return Int_Type is <>;
     with function "+"   (Left, Right : Int_Type) return Int_Type is <>;
  function Get_Hamming (Position : Positive) return Int_Type;

  function Get_Hamming (Position : Positive) return Int_Type is
     function Is_Hamming (Number : Int_Type) return Boolean is
        Temporary : Int_Type := Number;
     begin
        while Temporary mod Two = Zero loop
           Temporary := Temporary / Two;
        end loop;
        while Temporary mod Three = Zero loop
           Temporary := Temporary / Three;
        end loop;
        while Temporary mod Five = Zero loop
           Temporary := Temporary / Five;
        end loop;
        return Temporary = One;
     end Is_Hamming;
     Result   : Int_Type := One;
     Previous : Positive := 1;
  begin
     while Previous /= Position loop
        Result := Result + One;
        if Is_Hamming (Result) then
           Previous := Previous + 1;
        end if;
     end loop;
     return Result;
  end Get_Hamming;

  -- up to 2**32 - 1
  function Integer_Get_Hamming is new Get_Hamming
     (Int_Type => Integer,
      Zero     => 0,
      One      => 1,
      Two      => 2,
      Three    => 3,
      Five     => 5);

  -- up to 2**64 - 1
  function Long_Long_Integer_Get_Hamming is new Get_Hamming
     (Int_Type => Long_Long_Integer,
      Zero     => 0,
      One      => 1,
      Two      => 2,
      Three    => 3,
      Five     => 5);

begin

  Ada.Text_IO.Put ("1) First 20 Hamming numbers: ");
  for I in 1 .. 20 loop
     Ada.Text_IO.Put (Integer'Image (Integer_Get_Hamming (I)));
  end loop;
  Ada.Text_IO.New_Line;
  Ada.Text_IO.Put_Line ("2) 1_691st Hamming number: " &
     Integer'Image (Integer_Get_Hamming (1_691)));
  -- even Long_Long_Integer overflows here
  Ada.Text_IO.Put_Line ("3) 1_000_000st Hamming number: " &
     Long_Long_Integer'Image (Long_Long_Integer_Get_Hamming (1_000_000)));

end Hamming;</lang>

Output:

1) First 20 Hamming numbers:  1 2 3 4 5 6 8 9 10 12 15 16 18 20 24 25 27 30 32 36
2) 1_691 st Hamming number:  2125764000

Execution terminated by unhandled exception
Exception name: CONSTRAINT_ERROR
Message: hamming.adb:34 overflow check failed
Call stack traceback locations:
0x403212 0x402fd7 0x402a87 0x7f8b99517584 0x4026d7

ALGOL 68

Works with: Algol 68 Genie 1.19.0

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>

C#

<lang csharp>using System; using System.Collections; using System.Collections.Generic; using System.Linq;

public class LazyList : IEnumerable<int> {

   int First;
   Func<LazyList> Rest;

   public LazyList(int first, Func<LazyList> rest)
   {
       First = first;
       Rest = rest;
   }

   public LazyList Map(Func<int, int> function)
   {
       return new LazyList(function(First), () => Rest().Map(function));
   }

   public LazyList Merge(LazyList list)
   {
       if (this == null)
       {
           return list;
       }
       if (list == null)
       {
           return this;
       }
       if (First < list.First)
       {
           return new LazyList(First, () => Rest().Merge(list));
       }
       if (First > list.First)
       {
           return new LazyList(list.First, () => Merge(list.Rest()));
       }
       return new LazyList(First, () => Rest().Merge(list.Rest()));
   }

   public IEnumerator<int> GetEnumerator()
   {
       yield return First;
       if (Rest() != null)
       {
           foreach (var first in Rest())
           {
               yield return first;
           }
       }
   }

   IEnumerator IEnumerable.GetEnumerator()
   {
       return (IEnumerator)GetEnumerator();
   }

}

class Program {

   static void Main()
   {
       var hamming = default(LazyList);
       hamming = new LazyList(1, () => hamming.Map(x => x * 2).Merge(hamming.Map(x => x * 3).Merge(hamming.Map(x => x * 5))));
       foreach (var number in hamming.Take(20))
       {
           Console.WriteLine(number);
       }
   }

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

Translation of: Scala

<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 .

Translation of: Haskell

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 .

Fortran

Works with: Fortran version 90 and later

Using big_integer_module from here[1] <lang fortran>program Hamming_Test

 use big_integer_module
 implicit none
 
 call Hamming(1,20)
 write(*,*)
 call Hamming(1691)
 write(*,*)
 call Hamming(1000000)

contains

subroutine Hamming(first, last)

 integer, intent(in) :: first
 integer, intent(in), optional :: last
 integer :: i, n, i2, i3, i5, lim
 type(big_integer), allocatable :: hnums(:)

 if(present(last)) then
   lim = last
 else
   lim = first
 end if

 if(first < 1 .or. lim > 2500000 ) then
   write(*,*) "Invalid input"
   return
 end if
 
 allocate(hnums(lim))
 
 i2 = 1 ;  i3 = 1 ; i5 = 1  
 hnums(1) = 1
 n = 1
 do while(n < lim)
   n = n + 1
   hnums(n) = mini(2*hnums(i2), 3*hnums(i3), 5*hnums(i5))
   if(2*hnums(i2) == hnums(n)) i2 = i2 + 1
   if(3*hnums(i3) == hnums(n)) i3 = i3 + 1
   if(5*hnums(i5) == hnums(n)) i5 = i5 + 1
 end do
 
 if(present(last)) then
   do i = first, last
     call print_big(hnums(i))
     write(*, "(a)", advance="no") " "
   end do
 else
   call print_big(hnums(first))
 end if
 
 deallocate(hnums)

end subroutine

function mini(a, b, c)

 type(big_integer) :: mini
 type(big_integer), intent(in) :: a, b, c
  
 if(a < b ) then 
   if(a < c) then
     mini = a
   else
     mini = c
   end if
 else if(b < c) then
   mini = b
 else
   mini = c
 end if

end function mini end program</lang> Output

1 2 3 4 5 6 8 9 10 12 15 16 18 20 24 25 27 30 32 36 
2125764000
519312780448388736089589843750000000000000000000000000000000000000000000000000000000

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.

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

produce (the next) 3 elements by selection and multiplication:

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

<lang java> import java.math.BigInteger;

public class HammingRosetta { 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 h3[]=new BigInteger[limit/2]; final BigInteger h5[]=new BigInteger[limit];

BigInteger x = BigInteger.ONE; int j=0,k=0; int jj=0,kk=0;

BigInteger x2=TWO; BigInteger x3=THREE; BigInteger x5=FIVE;

for(int n=1; n<limit; n++) {

                       x=x2.min(x3.min(x5));

if (x2.equals(x)) { h3[jj++]=h5[kk++]=x; x2=TWO.multiply(x); }

else if (x3.equals(x)) { h3[jj++]=h5[kk++]=x; x3=THREE.multiply(h3[j]); h3[j++]=null; }

else { h5[kk++]=x; x5=FIVE.multiply(h5[k]); h5[k++]=null; } } return x; }

public static void main(String[] args) { System.out.println(hamming(Integer.valueOf(args[0]))); }

} </lang>

JavaScript

Works with: JavaScript version 1.7

Works with: Firefox version 2

Translation of: Ruby

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

Translation of: Haskell

<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, islice 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

print list(islice(raymonds_hamming(), 20)) print islice(raymonds_hamming(), 1689, 1690).next() print islice(raymonds_hamming(), 999999, 1000000).next()</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 $2^{31}$ . Works for larger limits, however to find the $1000000^{th}$ 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>

REXX

<lang rexx> numeric digits 100 call hamming 1,20 call hamming 1691 call hamming 1000000 exit

hamming: procedure; parse arg x,y if y== then y=x p2=1 p3=1 p5=1 a.=0 a.1=1

 do n=2 for y-1
 a.n=min(2*a.p2, 3*a.p3, 5*a.p5)
 if 2*a.p2==a.n then p2=p2+1
 if 3*a.p3==a.n then p3=p3+1
 if 5*a.p5==a.n then p5=p5+1
 end

say

 do j=x to y
 say 'hamming('j")=" a.j
 end

say say right('length of last hamming number='length(a.y),70) say return </lang> 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

                                       length of last hamming number=2


hamming(1691)= 2125764000

                                      length of last hamming number=10


hamming(1000000)= 519312780448388736089589843750000000000000000000000000000000000000000000000000000000

                                      length of last hamming number=84

Ruby

Translation of: Scala

This method

Works with: Ruby version 1.8.6

<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,

Works with: Ruby version 1.9.1

, 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

Works with: GNU 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.

Works with: Tcl version 8.6

<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 $n$ to the $n$ -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>