## 357 Numbers

### March 6, 2015

We first note that 1 is a member of the set. Although it is not divisible by 3, 5 or 7, neither is it divisible by any number other than 3, 5 or 7, so we will say that 1 is in the set. This conforms to the mathematical definition of the set { *x* = 3^{i} · 5^{j} · 7^{k} | *i*, *j*, *k* ≥ 0 }.

The straight forward solution is to pick a number and divide repeatedly by 3, 5 and 7; if the residue is 1, then the number is a member of the set:

`(define (357? n)`

(while (zero? (modulo n 3)) (set! n (/ n 3)))

(while (zero? (modulo n 5)) (set! n (/ n 5)))

(while (zero? (modulo n 7)) (set! n (/ n 7)))

(= n 1))

Then it is simple to check all the odd numbers less than the desired maximum number in the set:

`(define (list357 n) (filter 357? (range 1 n 2)))`

`> (time (length (list357 1000000)))`

(time (length (list357 1000000)))

29 collections

668 ms elapsed cpu time, including 54 ms collecting

675 ms elapsed real time, including 54 ms collecting

240014352 bytes allocated, including 245428480 bytes reclaimed

203

That works, but is rather slow, as it must examine half a million numbers and discard 499,797 of them. A better approach is to generate only those numbers that are in the set, which can be done by induction: 1 is in the set, and for any number *x* in the set, 3*x*, 5*x* and 7*x* are also in the set. So the function that computes the set pushes 1 into an initially-empty priority queue, then repeatedly extracts the smallest number in the priority queue, inserts three new numbers into the priority queue, and stops when requested:

`(define (list357 n)`

(let ((pq (pq-insert < 1 pq-empty)))

(let loop ((357s (list)))

(if (pq-empty? pq) (reverse 357s)

(let ((x (pq-first pq)))

(while (and (not (pq-empty? pq))

(= (pq-first pq) x))

(set! pq (pq-rest < pq)))

(when (< (* 3 x) n)

(set! pq (pq-insert < (* 3 x) pq)))

(when (< (* 5 x) n)

(set! pq (pq-insert < (* 5 x) pq)))

(when (< (* 7 x) n)

(set! pq (pq-insert < (* 7 x) pq)))

(loop (cons x 357s)))))))

That’s a lot faster:

`> (time (length (list357 1000000)))`

(time (length (list357 1000000)))

no collections

2 ms elapsed cpu time

2 ms elapsed real time

190704 bytes allocated

203

It can find the set to 10^{35} in the same time the brute-force method finds the set to a million:

`> (time (length (list357 #e1e35)))`

(time (length (list357 100000000000000000000000000000000000)))

6 collections

681 ms elapsed cpu time, including 18 ms collecting

685 ms elapsed real time, including 19 ms collecting

53436416 bytes allocated, including 48876240 bytes reclaimed

27608

And even finding the set to a googol isn’t too traumatic:

`> (time (length (list357 #e1e100)))`

(time (length (list357 10000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000)))

177 collections

19308 ms elapsed cpu time, including 842 ms collecting

19317 ms elapsed real time, including 843 ms collecting

1488283648 bytes allocated, including 1746698336 bytes reclaimed

609456

We implemented the priority queue using the pairing heaps of a previous exercise. You can run the program at http://ideone.com/YwnAQ8. You might also look at the algorithm of A108347, which is even worse than the slow method described above; it factors each number, forms the factors into a set, then compares the set to the set {3 5 7}.

But there must some conspiracy afoot, since my result is 203 :S

(length (sort (let ((max 1000000))

(loop

:for 3s = 1 then (* 3 3s)

:while (< 3s max)

:append (loop

:for 5s = 1 then (* 5 5s)

:for 35s = (* 3s 5s)

:while (< 35s max)

:append (loop

:for 7s = 1 then (* 7 7s)

:for 357s = (* 35s 7s)

:while (< 357s max)

:collect 357s))))

(function <)))

203

In Python. The first method (div357) is more or less brute force. It generates all possible products of powers of 3,5 and 7. The second method uses a queue to generate the same numbers, but this method is in Python a lot slower.

I keep getting 202 instead of 203, have gone over everything but can’t seem to find the 203rd one. Are you guys including 0 or 1 there by any chance? Can anyone publish a list of the numbers found?

@FelipeCeotto: The number 1 is included in the set. See the first paragraph of the solution page.

Thanks, @programmingpraxis

My somewhat verbose solution in C#:

Reminiscent of the Hamming Numbers exercise. My Haskell solution is similar to my submission for that exercise:

Python. Returns a sorted list of the numbers. Could be sped up by pulling powersof() out of the inner loops. As is, p357(1000000) only takes about a millisecond.

This is a Scheme implementation of Mike’s Python, with infrastructure. My generators are procedures that return successive values on successive calls and then raise an exception (or return the values from a sentinel thunk); their ability to return multiple values on each call is not used in this exercise. The mechanism was originally inspired by a post in the Python newsgroup. I’m not sure I’ve got it exactly right, but I’m also not aware of where any problems might lie :/

;;; (length (sorted (p357 1000000))) => 203

(import (rnrs exceptions (6))) ;; in Guile 2.0.5

(define (p357 limit)

(generator

(lambda (return)

(for (powers-of 3 limit)

(lambda (p3)

(for (powers-of 5 (ceiling-quotient limit p3))

(lambda (p5)

(for (powers-of 7 (ceiling-quotient limit (* p3 p5)))

(lambda (p7)

(return (* p3 p5 p7)))))))))))

(define (powers-of n limit)

(generator

(lambda (return)

(do ((p 1 (* p n)))

((>= p limit))

(return p)))))

(define (for gen proc)

(guard

(exn

((eq? exn 'the-end-of-generator)))

(call-with-values gen proc)

(for gen proc)))

(define (sorted gen)

(let ((result '()))

(guard

(exn

((eq? exn 'the-end-of-generator) (sort result <)))

(let loop ()

(set! result (cons (gen) result))

(loop)))))

`(define generator`

(case-lambda

((w)

(let ((s (lambda ()

(raise 'the-end-of-generator))))

(generator w s)))

((w s)

(let* ((receiver values)

(terminator values)

(condition '())

(body (lambda (yield)

(w yield)

(call-with-values s receiver)))

(yield (lambda vs

(call-with-current-continuation

(lambda (k)

(set! body k)

(apply receiver vs))))))

(lambda ()

(call-with-current-continuation

(lambda (k)

(set! receiver k)

(call-with-current-continuation

(lambda (k)

(set! terminator k)

(guard

(con (else

(set! condition con)

(terminator)))

(body yield))

(terminator)))

(raise condition))))))))

@Jussi: You could also use the generators defined in the Standard Prelude:

`> (define-generator (make-gen357)`

(let ((pq (pq-insert < 1 pq-empty)))

(let loop ()

(let ((x (pq-first pq)))

(while (and (not (pq-empty? pq))

(= (pq-first pq) x))

(set! pq (pq-rest < pq)))

(set! pq (pq-insert < (* 3 x) pq))

(set! pq (pq-insert < (* 5 x) pq))

(set! pq (pq-insert < (* 7 x) pq))

(yield x) (loop)))))

> (define gen357 (make-gen357))

> (gen357)

1

> (gen357)

3

> (gen357)

5

> (gen357)

7

> (gen357)

9

> (gen357)

15

> (gen357)

21

> (define gen357-2 (make-gen357))

> (gen357-2)

1

> (gen357-2)

3

> (gen357-2)

5

> (gen357-2)

7

> (gen357-2)

9

> (gen357)

25

> (gen357)

27

> (gen357)

35

@Praxis, I can’t get your Standard Prelude define-generator to do what I want. There seems to be some problem defining generators in terms of generators, which is precisely what I wanted to emulate Mike’s code. This is again in Guile 2.0.5, with your definition from Standard Prelude, and my power generator modified to raise the specific condition when it runs out of values, as follows.

(import (rnrs exceptions (6)))

;;; ...

(define-generator (powers-of n limit)

(do ((p 1 (* p n)))

((>= p limit) (raise 'the-end-of-generator))

(yield p)))

Now, this works (with my little looping construct as before):

(define (e357 limit)

(for (powers-of 3 limit)

(lambda (p3)

(for (powers-of 5 (ceiling-quotient limit p3))

(lambda (p5)

(for (powers-of 7 (ceiling-quotient limit (* p3 p5)))

(lambda (p7)

(write (* p3 p5 p7))

(write-char #\space))))))))

That’s an “ordinary” procedure that simply writes out the products from the value-generating procedures. For example:

scheme@(guile-user)> (e357 10)

1 7 5 3 9 $1 = #t

The final part is Guile telling that the expression returned a truth value. But

I didn’t want to write those products out. I wanted to generate them to be sorted and counted. And this fails:

(define-generator (p357 limit)

(for (powers-of 3 limit)

(lambda (p3)

(for (powers-of 5 (ceiling-quotient limit p3))

(lambda (p5)

(for (powers-of 7 (ceiling-quotient limit (* p3 p5)))

(lambda (p7)

(yield (* p3 p5 p7)))))))))

It yields the first product all right, but then it fails to continue:

scheme@(guile-user)> (let ((g (p357 10))) (write (g)) (write (g)))

1ERROR: In procedure #<continuation 92456f0>:

ERROR: Too few values returned to continuation

I’ve reduced the failure down to what I thought would be two identity mappings:

(define-generator (gd-id gen) (let loop () (yield (gen)) (loop)))

(define-generator (bd-id gen) (for gen (lambda (v) (yield v))))

The first works. The second fails after the first value:

scheme@(guile-user)> (for (gd-id (powers-of 3 10)) write)

139$1 = #t

scheme@(guile-user)> (for (bd-id (powers-of 3 10)) write)

1ERROR: In procedure #<continuation 9c79220>:

ERROR: Too few values returned to continuation

It’s the second one that is like the nested loops in my Mike-emulator.

In my own system the following three all appear to work the same, as I expected.

(define (gd-id gen) (generator (lambda (yield) (let loop () (yield (gen)) (loop)))))

(define (od-id gen) (generator (lambda (yield) (for gen yield))))

(define (wd-id gen) (generator (lambda (yield) (for gen (lambda (v) (yield v))))))

So I failed to use the Standard Prelude in place of my own generators. Did I try to do something that’s wrong? Much apologies if I’m making some very stupid mistake!

Inspired by Graham’s Haskell code, I came up with these python generators. s357() generates the values in order without repetition.

@Mike’s Python code above is quite fast (for Python) due to the efficiency of the heapq merge function, but Python doesn’t really like recursion very much as it isn’t optimized for tail recursion as can be shown by changing the definition of just one of the series generators, the s3() one, to a version not requiring iterator recursion which speeds it up by almost a factor of two, as follows:

The recursion is slower due to the overhead of successive iterator function calls to implement it that way.

Now we could re-write the algorithm to not require the recursive iteration, but as Python doesn’t optimize successive calls very well or inline code, there are even more nested function calls and the code is slightly slower. The following code uses Co-Inductive Streams implemented as two-tuples with the first element carrying the value and the second a “thunk’ lambda for the continuation function to the next value:of the “stream”:

For this algorithm, it is better to use a language optimized for tail calls such as Scheme, as follows (this implementation uses Co-Inductive Streams implemented as a Pair, with the elements as per the two-tuple of the Python code):

(define (s357)

(define (merge x y)

(let ((xv (car x)) (yv (car y)))

(cond (( xv yv) (cons yv (lambda () (merge x ((cdr y))))))

(else (cons xv (lambda () (merge ((cdr x)) y)))))))

(define (smult s m) (cons (* (car s) m) (lambda () (smult ((cdr s)) m))))

(define (s3 n) (cons n (lambda () (s3 (* n 3)))))

(define (s35) (cons 1 (lambda () (merge (s3 3) (smult (s35) 5)))))

(let ((s35r (cdr (s35)))) (cons 1 (lambda () (merge (s35r) (smult (s357) 7))))))

;;; testing…

(define limit 100000000000000000000000000000000000)

(define (countemto lmt)(do ((nxt (s357) ((cdr nxt))) (cnt 0 (+ cnt 1))) ((> (car nxt) lmt) cnt)))

(display (time (countemto limit))) (newline)

where even run as interpreted code on http://codepad.org/DUIJHxBA, it runs faster than Python as shown following:

cpu time: 890 real time: 965 gc time: 501

27608

It will run 10’s of times faster compiled to machine code with Gambit-C, Bigloo, or Chez Scheme compilers with optimizations and unsafe modes turned on.

Note that much of the time used for the codepad interpreter is in garbage collecting (the codepad collector must be pretty poor); for the 64-bit Petite Chez Scheme, the output is as follows (note that the main difference is the amount of time spent in garbage collecting, with Chez using very little):

The Scheme listing didn’t come out correctly the first time, although it is correct at the codepage link: http://codepad.org/DUIJHxBA.

The code listing should be as follows:

Another merging scheme is possible, so there are no duplicates, as discussed at http://stackoverflow.com/q/12480291/849891. In Haskell:

h357 = 1 : q3 — length . takeWhile ( 203

where

q7 = 7 : map (* 7) q7

q5 = 5 : merge q7 (map (* 5) q5)

q3 = 3 : merge q5 (map (* 3) q3)

It can also be written as a one-liner, as

1:foldr (\n s->fix ((n:) . merge s . map (n*))) [] [3,5,7]

Due to Rosettacode’s user “Ledrug”.

Sorry for messing up the markup. I should have also read through to the end, a similar scheme was already proposed above.

It can also be written as a one-liner, as

Due to Rosettacode’s user “Ledrug”.

Great contribution by @Will_Ness as to the significance of these solutions!

As per whatever @Mike based the Python solution on, we were working toward the Daniel Fischer solution from the stackoverflow link.

The “one-liner” (well, not quite as it still requires the “merge” clause) Haskell solution isn’t possible in Scheme (or Python) as there is no built-in lazy lists and therefore no fixpoint recursion; moreover, reversing the order from 3/5/7 to 7/5/3 doesn’t have any benefit and in fact makes the current Scheme (and Python) code slower for larger ranges due to not taking as much advantage of the simplified sequence for the most often repeated ‘3’ series: Also, using a direct sequence multiply function instead of a map applied to a (lazy) list is more efficient being more direct. Finally, @Mike’s use of Python iterator sequences or my use of Co-Inductive Streams (CIS’s) rather than some sort of Haskell-like lazy lists (Streams as per Schemes SRFI-41, authored by @praxis) limits the speed of our current implementations due to having to recalculate previous values of the sequence due to no memoization being applied; however, the speed of the operations and less demands on garbage collection tends to cancel the theoretical advantage out for these languages.

In order to try to take full advantage of the Daniel Fischer discovery and the constant factor advantage, a memoized lazy stream is needed, which is easy to implement from the CIS code by just changing the lambda thunks to use the R5rS delay/force primitives (for this limited purpose, there is no need to use the full SRFI-18 implementation). The code which is more or less equivalent to @Will_Ness’s Haskell code then is as follows:

which produces the best time with the Gambit-C compiler using optimized and unsafe code:

The above code is slower than the previous code for most Scheme implementations and slower than compiled Haskell due to not optimizing the lazy stream operations and the greater demands on (likely less efficient) garbage collection; as well, the reversing of the prime processing actually takes longer due to increased merge operations for the above reasons.

Live and learn!

This is my implementation in Java.

Is this what you mean?

public class DivisibleBy357Only

{

private static final int MAX_NUM = 1000000;

public static Set findNumDivisibleBy357Only()

{

Set divisibleNum = new HashSet();

Queue ququeNum = new LinkedList();

ququeNum.add(1);

Integer curQueueNum, product;

int[] baseNum = {3, 5, 7};

while( (curQueueNum = ququeNum.poll()) != null)

{

for(int i=0; i < baseNum.length; i++)

{

product = curQueueNum * baseNum[i];

if( product < MAX_NUM )

{

ququeNum.add(product);

divisibleNum.add(product);

}

}

}

return divisibleNum;

}

}

@Xin: Your Java program does not produce the 357factors in order, meaning if that were required you would have to sort it using a TreeSet rather than a HashSet at a cost in performance for larger ranges (O(log n) per entry rather than O(1))…

I couldn’t resist a C++ version of the lazy list merging solution with 3 lists, using STL deques:

Perhaps not as elegant as the Haskell solution and certainly not a one liner. Perf tells me it take 1.8M cycles to count the 203 numbers less that 1000000.

Here is a recursive solution in c++. May stack frame depth will only be log 3 1000000 = 12 stack frames.

@brooknovak: Nice solution – we don’t need to fully recalculate the sum each time though (and people might worry about floating point equality) and we don’t need to recurse for numbers we have seen before, so (with C++11 for auto) we can just do:

And a variant on the generating all triples approach, with no nested loops. Numbers not generated in order but that doesn’t matter if we are just counting them:

It just occurred to me that t3 serves no purpose in that code and can be removed.

And here is a generalization:

I keep coming back to this exercise wondering why Haskell is tens of times faster than any of the other functional languages and algorithms used here so did some investigating and have come up with the following:

1. It has nothing to do with its conciseness (although it definitely is that) as it is just as fast when the algorithm is written in the same form as Scheme (as I showed in a couple of examples in the comments above) or in F# (not presented here – using sequences or better Co Inductive Streams, as the results were the same).

2. I ran the examples under Haskell as @Will_Ness represented them (as he obtained them from his sources) and they are indeed about 60 times faster than as compiled under the best implementation of Scheme that generates efficient Cee code that I have at my disposal (Bigloo or Gambit-C) or F#, which is about the same speed as compiled Scheme. I think that the Haskell code (using the LLVM back end) is faster for two reasons, as follows:

a) Haskell uses that it is a pure hygienic functional language with generally nothing but immutable state to tune its custom garbage collector for very high efficiency as compared to the (Boehm) garbage collector used by Bigloo or even the custom garbage collector used by Gambit-C since Scheme does not enforce purity as Haskell does (most programmer use this in lots of “set!” procedures to modify values and variables), and

b) the GHC Haskell compiler (to some extent in combination with the LLVM optimization process) is extremely efficient at optimizing the output machine code as to eliminating procedure/function calls and converting them to simple loops and inline code. The Scheme compilers (and the F# compilation process) do some optimization along these lines in determining code that can be inlined and in tail call optimizations as to producing loops, but aren’t very efficient when it comes to iteration (Scheme doesn’t have built in iterators other than lists, which are not lazy and therefore not true iterators) and F# has lazy sequences which are very poorly optimized, using many many nested procedure calls to implement. In particular, Haskell’s lazy lists are implemented in such a way that the compiler can turn the merge function into a simple comparison and a change of link destination rather than a procedure call to a “thunk” (zero argument function) which in turn calls another function with its closure parameters, at 10’s of CPU clock cycles per function call.

In short, it appears Haskell can’t be beat in performance nor in conciseness for this sort of reasonably complex recursive list processing due to its pure functional nature, and simplistic languages/compilers such as F#/Scheme can’t touch it unless their compilers invested the years of research required in garbage collection optimizations, strictness analysis, and pure functional optimizations.

BTW, while Daniel Fischer’s optimizations to reduce redundant work does cause the large constant factor in better performance for the algorithm no matter the implementation language, changing the order of operations doesn’t have that large an effect for these 3/5/7 combination (perhaps it has more for the 2/3/5 combination he tested), and @Will_Ness’s “one-liner” is actually slightly faster when written as follows (written as a function so it will reevaluate when run under the REPL):

The GHC Haskell compiler is actually generating similar machine code as do the latest C++11 codes from @matthew except that it doesn’t use mutable state variables and thus there is considerable garbage collection activity taking slightly longer than the actual computation time; that said, the GHC Haskell garbage collector is extremely optimized for there being no mutable state and thus is more efficient than garbage collectors/automatic reference counting that must deal with that so is faster than C++11’s automatic safe_pointer handling if that were used.

There is also that the multiple precision (bitnum) numeric library Haskell uses is written in Haskell and thus much faster than the (Cee) GMP library used by most Scheme imlementations. The C++ code others have supplied for fixed precision is likely faster, but may well be slower if written for big numbers, as in general the routines do not eliminate repeating occurrences such as 3 time 5 and 5 times 3 also occurring in many combinations of factors, whereas the @Will_Ness supplied references and code as well as the @Mike Python code and my Scheme codes eliminate redundancies; for this reason the Haskell code may well be faster than Cee in spite of using immutable state and thus requiring the garbage collections.

Interesting analysis, thanks. Some impressive optimizations by the Haskell compiler.

@matthew, yes, Haskell is fast and often as fast or faster than C, especially when using the LLVM backend. In this case, the algoirthm is well suited to lazy list processing, which in the merge operation removes redundant calculaatons and leaves the output sorting in ascending order for free in the same operations. However, the whole think also requires the memoization of lazy list processing, which I missed in my earlier posts; the problem with those (Scheme) posts is that the internal lists are used as procedures, which means that when they are used recursively as required by the last two streams, they require that all calculations must be repeated, meaning that instead of one comparison and just over one multiplication per output number, there are about forty times that many.

The code uses a stream/lazy list implementation based on pairs and the R6RS delay/force promise handling procedures, with internal recursive procedures changes to stream variables so that they don’t continuously get re-calculated.

The following code finds the first 10,000,000 (ten million) numbers in the sequence in about 7.6 seconds including about 2.1 seconds for garbage collection using compiled Gambit-C with maximum optimizations turned on, which is only a little over twice as slow as the best optimized Haskell code:

The above code interpreted under Petite Chez Scheme (as praxis uses) calculates the sequence for a million numbers in about 2.8 seconds with about 2.0 seconds used for garbage collection and the sequence up to a googol in about 1.8 seconds with about 1.0 seconds used for garbage collection.

It calculates the first one million numbers in the sequence in about 2.2 seconds with about 0.5 seconds used for garbage collection under codepad: http://codepad.org/W3QrjxmS.

With DrRacket using R6RS format and with all debugging turned off, it calculates the first million numbers in about 3.6 seconds with about 2.4 seconds used for garbage collection.

The Haskell code posted earlier is still faster than the fastest Scheme implementation of the above code due to its rule-based optimization specializations for lazy list by a factor of over two, but if the Haskell code is written to emulate the Scheme code using roll-your-own streams with thunks used as a delayed calculation means so that the specialized optimizations aren’t used, then the code is about the same speed as the very fastest Scheme implementation (Gambit-C).

One needs to remember that memorization if often desired to avoid repeated recursive calculations, such as required here for the second two sequences that refer to their own previous results. The recursive references to earlier values in the list mean that memory residency will be increased, but the huge increase in improved execution speed is likely worth it.