Iterated digits squaring
If you add the square of the digits of a Natural number (an integer bigger than zero), you always end with either 1 or 89:
15 -> 26 -> 40 -> 16 -> 37 -> 58 -> 89 7 -> 49 -> 97 -> 130 -> 10 -> 1
An example in Python:
<lang python>>>> step = lambda x: sum(int(d) ** 2 for d in str(x)) >>> iterate = lambda x: x if x in [1, 89] else iterate(step(x)) >>> [iterate(x) for x in xrange(1, 20)] [1, 89, 89, 89, 89, 89, 1, 89, 89, 1, 89, 89, 1, 89, 89, 89, 89, 89, 1]</lang>
- Task
- Count how many number chains for integers 1 <= n < 100_000_000 end with a value 89.
Or, for much less credit - (showing that your algorithm and/or language is slow):
- Count how many number chains for integers 1 <= n < 1_000_000 end with a value 89.
This problem derives from the Project Euler problem 92.
- Cf
D
A simple memoizing partially-imperative brute-force solution: <lang d>import std.stdio, std.algorithm, std.range, std.functional;
uint step(uint x) pure nothrow @safe @nogc {
uint total = 0;
while (x) {
total += (x % 10) ^^ 2;
x /= 10;
}
return total;
}
uint iterate(in uint x) nothrow @safe {
return (x == 89 || x == 1) ? x : x.step.memoize!iterate;
}
void main() {
iota(1, 100_000_000).filter!(x => x.iterate == 89).count.writeln;
}</lang>
- Output:
85744333
The run-time is about 10.9 seconds compiled with ldc2.
A fast imperative brute-force solution: <lang d>void main() nothrow @nogc {
import core.stdc.stdio: printf;
enum uint magic = 89; enum uint limit = 100_000_000; uint[(9 ^^ 2) * 8 + 1] lookup = void;
uint[10] squares;
foreach (immutable i, ref x; squares)
x = i ^^ 2;
foreach (immutable uint i; 1 .. lookup.length) {
uint x = i;
while (x != magic && x != 1) {
uint total = 0;
while (x) {
total += squares[(x % 10)];
x /= 10;
}
x = total;
}
lookup[i] = x == magic; }
uint magicCount = 0;
foreach (immutable uint i; 1 .. limit) {
uint x = i;
uint total = 0;
while (x) {
total += squares[(x % 10)];
x /= 10;
}
magicCount += lookup[total]; }
printf("%u\n", magicCount);
}</lang> The output is the same. The run-time is less than 3 seconds compiled with ldc2.
A more efficient solution: <lang d>import core.stdc.stdio, std.algorithm, std.range;
enum factorial = (in uint n) pure nothrow @safe @nogc
=> reduce!q{a * b}(1u, iota(1u, n + 1));
uint iLog10(in uint x) pure nothrow @safe @nogc in {
assert(x > 0);
} body {
return (x >= 1_000_000_000) ? 9 :
(x >= 100_000_000) ? 8 :
(x >= 10_000_000) ? 7 :
(x >= 1_000_000) ? 6 :
(x >= 100_000) ? 5 :
(x >= 10_000) ? 4 :
(x >= 1_000) ? 3 :
(x >= 100) ? 2 :
(x >= 10) ? 1 : 0;
}
uint nextStep(uint x) pure nothrow @safe @nogc {
typeof(return) result = 0;
while (x > 0) {
result += (x % 10) ^^ 2;
x /= 10;
}
return result;
}
uint check(in uint[] number) pure nothrow @safe @nogc {
uint candidate = reduce!((tot, n) => tot * 10 + n)(0, number);
while (candidate != 89 && candidate != 1)
candidate = candidate.nextStep;
if (candidate == 89) {
uint[10] digitsCount;
foreach (immutable d; number)
digitsCount[d]++;
return reduce!((r, c) => r / c.factorial)
(number.length.factorial, digitsCount);
}
return 0;
}
void main() nothrow @nogc {
enum uint limit = 100_000_000; immutable uint cacheSize = limit.iLog10;
uint[cacheSize] number; uint result = 0; uint i = cacheSize - 1;
while (true) {
if (i == 0 && number[i] == 9)
break;
if (i == cacheSize - 1 && number[i] < 9) {
number[i]++;
result += number.check;
} else if (number[i] == 9) {
i--;
} else {
number[i]++;
number[i + 1 .. $] = number[i];
i = cacheSize - 1;
result += number.check;
}
}
printf("%u\n", result);
}</lang> The output is the same. The run-time is about 0.04 seconds or less. This third version was ported to D and improved from: mathblog.dk/project-euler-92-square-digits-number-chain/
Java
<lang java>import java.util.stream.IntStream;
public class IteratedDigitsSquaring {
public static void main(String[] args) {
long r = IntStream.range(1, 100_000_000)
.parallel()
.filter(n -> calc(n) == 89)
.count();
System.out.println(r);
}
private static int calc(int n) {
while (n != 89 && n != 1) {
int total = 0;
while (n > 0) {
total += Math.pow(n % 10, 2);
n /= 10;
}
n = total;
}
return n;
}
}</lang>
85744333
Perl
<lang perl>use warnings; use strict;
my @sq = map { $_ ** 2 } 0 .. 9; my %cache; my $cnt = 0;
sub Euler92 {
my $n = 0 + join( , sort split( , shift ) );
$cache{$n} //= ($n == 1 || $n == 89) ? $n :
Euler92( sum( @sq[ split , $n ] ) )
}
sub sum {
my $sum; $sum += shift while @_; $sum;
}
for (1 .. 100_000_000) {
++$cnt if Euler92( $_ ) == 89;
}
print $cnt;</lang>
85744333
Perl 6
This fairly abstract version does caching and filtering to reduce the number of values it needs to check and moves calculations out of the hot loop, but is still interminably slow... even for just up to 1,000,000.
<lang perl6>constant @sq = ^10 X** 2; my $cnt = 0; my %cache;
sub Euler92($n) {
%cache{$n} //=
$n == any(1,89) ?? $n !!
Euler92( [+] @sq[$n.comb] )
}
for 1 .. 1_000_000 -> $n {
my $i = +$n.comb.sort.join; ++$cnt if Euler92($i) == 89;
}
say $cnt;</lang>
- Output:
856929
All is not lost, however. Through the use of gradual typing, Perl 6 scales down as well as up, so this jit-friendly version is performant enough to brute force the larger calculation: <lang perl6>my @cache; @cache[1] = 1; @cache[89] = 89;
sub Euler92(int $n) {
$n < 649 # 99,999,999 sums to 648, so no point remembering more
?? (@cache.at_pos($n) //= ids($n))
!! ids($n)
}
sub ids(int $num --> int) {
my int $n = $num;
my int $ten = 10;
my int $sum = 0;
my int $t;
my int $c;
repeat until $n == 89 or $n == 1 {
$sum = 0;
repeat {
$t = $n div $ten;
$c = $n - $t * $ten;
$sum = $sum + $c * $c;
} while $n = $t;
$n = @cache.at_pos($sum) // $sum;
}
$n;
}
my int $cnt = 0; for 1 .. 100_000_000 -> int $n {
$cnt = $cnt + 1 if Euler92($n) == 89;
} say $cnt;</lang>
- Output:
85744333
This runs in under ten minutes. We can reduce this to 4 minutes by writing in the NQP (Not Quite Perl6) subset of the language: <lang perl6>my $cache := nqp::list_i(); nqp::bindpos_i($cache, 650, 0); nqp::bindpos_i($cache, 1, 1); nqp::bindpos_i($cache, 89, 89);
sub Euler92(int $n) {
$n < 650
?? nqp::bindpos_i($cache,$n,ids($n)) !! ids($n) }
sub ids(int $num --> int) {
my int $n = $num;
my int $ten = 10;
my int $sum = 0;
my int $t;
my int $c;
repeat until $n == 89 or $n == 1 {
$sum = 0; repeat { $t = nqp::div_i($n, $ten); $c = $n - $t * $ten; $sum = $sum + $c * $c; } while $n = $t; $n = nqp::atpos_i($cache,$sum) || $sum;
} $n;
}
my int $cnt = 0; for 1 .. 100_000_000 -> int $n {
$cnt = $cnt + 1 if Euler92($n) == 89;
} say $cnt;</lang>
Python
Python: Combinatorics
<lang python>from math import ceil, log10, factorial
def next_step(x):
result = 0
while x > 0:
result += (x % 10) ** 2
x /= 10
return result
def check(number):
candidate = 0
for n in number:
candidate = candidate * 10 + n
while candidate != 89 and candidate != 1:
candidate = next_step(candidate)
if candidate == 89:
digits_count = [0] * 10
for d in number:
digits_count[d] += 1
result = factorial(len(number))
for c in digits_count:
result /= factorial(c)
return result
return 0
def main():
limit = 100000000 cache_size = int(ceil(log10(limit))) assert 10 ** cache_size == limit
number = [0] * cache_size result = 0 i = cache_size - 1
while True:
if i == 0 and number[i] == 9:
break
if i == cache_size - 1 and number[i] < 9:
number[i] += 1
result += check(number)
elif number[i] == 9:
i -= 1
else:
number[i] += 1
for j in xrange(i + 1, cache_size):
number[j] = number[i]
i = cache_size - 1
result += check(number)
print result
main()</lang>
- Output:
85744333
The run-time is less than half a second.
Python: Simple caching
<lang python>>>> from functools import lru_cache >>> @lru_cache(maxsize=1024) def ids(n): if n in {1, 89}: return n else: return ids(sum(int(d) ** 2 for d in str(n)))
>>> ids(15)
89
>>> [ids(x) for x in range(1, 21)]
[1, 89, 89, 89, 89, 89, 1, 89, 89, 1, 89, 89, 1, 89, 89, 89, 89, 89, 1, 89]
>>> sum(ids(x) == 89 for x in range(1, 100000000))
85744333
>>> </lang>
This took a much longer time, in the order of hours.
Python: Enhanced caching
Notes that the order of digits in a number does not affect the final result so caches the digits of the number in sorted order for more hits. <lang python>>>> from functools import lru_cache >>> @lru_cache(maxsize=1024) def _ids(nt): if nt in {('1',), ('8', '9')}: return nt else: return _ids(tuple(sorted(str(sum(int(d) ** 2 for d in nt)))))
>>> def ids(n):
return int(.join(_ids(tuple(sorted(str(n))))))
>>> ids(1), ids(15) (1, 89) >>> [ids(x) for x in range(1, 21)] [1, 89, 89, 89, 89, 89, 1, 89, 89, 1, 89, 89, 1, 89, 89, 89, 89, 89, 1, 89] >>> sum(ids(x) == 89 for x in range(1, 100000000)) 85744333 >>> _ids.cache_info() CacheInfo(hits=99991418, misses=5867462, maxsize=1024, currsize=1024) >>> </lang>
This took tens of minutes to run.
Ruby
<lang ruby>
- Count how many number chains for integers 1 < n < 100_000_000 end with a value 89.
- Nigel_Galloway
- August 26th., 2014.
require 'benchmark' F = Array.new(9){|n| n==0?1:(1..n).inject(:*)} N = (G=Array.new(649){|n| (i=((n%1000)/100)**2+((n%100)/10)**2+(n%10)**2)==89 ? 0:i}).collect{|n| while n>1 do n = G[n] end; n }
z = [0,0] t = Benchmark.measure do (0..9).each{|n| (n..9).each{|i| (i..9).each{|g| (g..9).each{|e| (e..9).each{|l| (l..9).each{|a| (a..9).each{|o| (o..9).each{|w|
next if w == 0
gn = 40320
nn = [0,0,0,0,0,0,0,0,0,0]
ng = [n,i,g,e,l,a,o,w]
ng.each{|n| nn[n] += 1}
nn.each{|n| gn = gn/F[n]}
z[N[ng.collect{|n| n*n}.inject(:+)]] += gn
}}}}}}}} end
puts "\n\n#{z[1]} numbers produce 1 and #{z[0]} numbers produce 89"
puts "\n\nTiming\n#{t}" </lang> {out}
14255666 numbers produce 1 and 85744333 numbers produce 89 Timing 0.250000 0.000000 0.250000 ( 0.245309)
zkl
Using brute force is a never ending process so need to be clever, which takes under a second.
<lang zkl>fcn check(number){ // a list of digits: 13 is L(0,0,0,0,0,0,1,3)
candidate:=number.reduce(fcn(sum,n){ sum*10 + n },0); // digits to int
while(candidate != 89 and candidate != 1) // repeatedly sum squares of digits
{ candidate = candidate.toString().reduce(fcn(sum,c){ sum + 1*c*c },0); }
if(candidate == 89){ // count permutations of these digits, they all sum to 89
digitsCount:=List(0,0,0,0,0,0,0,0,0,0);
foreach d in (number){ digitsCount[d] += 1; }
return(digitsCount.reduce(fcn(r,c){ r/factorial(c) },cacheBang)); // cacheBang==number.len()!
}
0 // this number doesn't sum to 89 (ie sums to 1)
} fcn factorial(n) { (1).reduce(n,fcn(N,n){ N*n },1) }
limit:=0d100_000_000; cacheSize:=limit.toFloat().log10().ceil().toInt(); number:=(0).pump(cacheSize,List().write,0); // list of zeros result:=0; i:=cacheSize - 1; var cacheBang=factorial(cacheSize); //== number.len()!
while(True){ // create numbers s.t. no set of digits is repeated
if(i == 0 and number[i] == 9) break;
if(i == cacheSize - 1 and number[i] < 9){ number[i] += 1; result += check(number); }
else if(number[i] == 9) i -= 1;
else{
number[i] += 1;
foreach j in ([i + 1 .. cacheSize - 1]){ number[j] = number[i]; }
i = cacheSize - 1;
result += check(number);
}
} println(result);</lang>
- Output:
85744333