Self numbers
You are encouraged to solve this task according to the task description, using any language you may know.
A number n is a self number if there is no number g such that g + the sum of g's digits = n. So 18 is not a self number because 9+9=18, 43 is not a self number because 35+5+3=43.
The task is:
Display the first 50 self numbers; I believe that the 100000000th self number is 1022727208. You should either confirm or dispute my conjecture.
224036583-1 is a Mersenne prime, claimed to also be a self number. Extra credit to anyone proving it.
- See also
AppleScript
I couldn't follow the math in the Wikipedia entry, nor the discussion and code here so far. But an initial expedient of generating a list of all the integers from 1 to just over ten times the required number of results and then deleting those that could be derived using the described method revealed the sequencing pattern on which the code below is based. The results match those from the "expedience" code as far as that can sensibly be pushed and the hundred millionth self number is reported as the hoped-for figure, so the intervening results are presumed to be correct too. The code would need modification (and a faster computer!) if self numbers with more than about 11 or 12 digits were required.
<lang applescript>(*
Base-10 self numbers by index (single or range). After the initial single-digit numbers, the sequence follows an observed pattern whereby self numbers occur in groups whose consecutive members are 11 apart. The interval between these groups is 2. They occur in blocks of ten: eight ten-number groups sandwiched between two shorter groups. The intervals between blocks and the lengths of the shorter groups depend the most signficant digit to have changed in the number sequence since they were last set.
- )
on selfNumbers(indexRange)
-- Subscript to monitor progress, store required results, and indicate when finished. set indexRange to indexRange as list script storer property startIndex : beginning of indexRange property endIndex : end of indexRange property counter : 0 property output : {} on doneAfter(thisSelf) set counter to counter + 1 if (counter < startIndex) then return false set end of my output to thisSelf return (counter = endIndex) end doneAfter end script -- Start with the single-digit odds. set thisSelf to -1 repeat 4 times set thisSelf to thisSelf + 2 if (storer's doneAfter(thisSelf)) then return storer's output end repeat -- Main loop. set shortGroupLength to 10 repeat -- Per block. -- First short group. Its length in the first block is 10, thereafter the same as that of the second -- short group of the preceding block. The numeric interval preceding it is related to this length. set thisSelf to thisSelf + 2 + (10 - shortGroupLength) * 13 if (storer's doneAfter(thisSelf)) then return storer's output -- Consecutive numbers within groups are 11 apart. repeat (shortGroupLength - 1) times set thisSelf to thisSelf + 11 if (storer's doneAfter(thisSelf)) then return storer's output end repeat -- Eight ten-number groups, each preceded by an interval of 2. repeat 8 times set thisSelf to thisSelf + 2 if (storer's doneAfter(thisSelf)) then return storer's output repeat 9 times set thisSelf to thisSelf + 11 if (storer's doneAfter(thisSelf)) then return storer's output end repeat end repeat -- Second short group. Its length depends on the most significant digit to tick over -- at the thousands boundary within it. The interval preceding it is 2. set thisSelf to thisSelf + 2 if (storer's doneAfter(thisSelf)) then return storer's output set i to 2 set shortGroupLength to 9 repeat until (i > shortGroupLength) set thisSelf to thisSelf + 11 if (storer's doneAfter(thisSelf)) then return storer's output if (thisSelf mod 1000 < 11) then set temp to thisSelf div 1000 repeat while (temp mod 10 is 0) set shortGroupLength to shortGroupLength - 1 set temp to temp div 10 end repeat end if set i to i + 1 end repeat end repeat return storer's output
end selfNumbers
-- Demo calls: -- First to fiftieth self numbers. selfNumbers({1, 50}) --> {1, 3, 5, 7, 9, 20, 31, 42, 53, 64, 75, 86, 97, 108, 110, 121, 132, 143, 154, 165, 176, 187, 198, 209, 211, 222, 233, 244, 255, 266, 277, 288, 299, 310, 312, 323, 334, 345, 356, 367, 378, 389, 400, 411, 413, 424, 435, 446, 457, 468}
-- One hundred millionth: selfNumbers(100000000) --> {1022727208}</lang>
C
Sieve based
About 25% faster than Go (using GCC 7.5.0 -O3) mainly due to being able to iterate through the sieve using a pointer. <lang c>#include <stdio.h>
- include <stdlib.h>
- include <time.h>
typedef unsigned char bool;
- define TRUE 1
- define FALSE 0
- define MILLION 1000000
- define BILLION 1000 * MILLION
- define MAX_COUNT 2*BILLION + 9*9 + 1
void sieve(bool *sv) {
int n = 0, s[8], a, b, c, d, e, f, g, h, i, j; for (a = 0; a < 2; ++a) { for (b = 0; b < 10; ++b) { s[0] = a + b; for (c = 0; c < 10; ++c) { s[1] = s[0] + c; for (d = 0; d < 10; ++d) { s[2] = s[1] + d; for (e = 0; e < 10; ++e) { s[3] = s[2] + e; for (f = 0; f < 10; ++f) { s[4] = s[3] + f; for (g = 0; g < 10; ++g) { s[5] = s[4] + g; for (h = 0; h < 10; ++h) { s[6] = s[5] + h; for (i = 0; i < 10; ++i) { s[7] = s[6] + i; for (j = 0; j < 10; ++j) { sv[s[7] + j+ n++] = TRUE; } } } } } } } } } }
}
int main() {
int count = 0; clock_t begin = clock(); bool *p, *sv = (bool*) calloc(MAX_COUNT, sizeof(bool)); sieve(sv); printf("The first 50 self numbers are:\n"); for (p = sv; p < sv + MAX_COUNT; ++p) { if (!*p) { if (++count <= 50) printf("%ld ", p-sv); if (count == 100 * MILLION) { printf("\n\nThe 100 millionth self number is %ld\n", p-sv); break; } } } free(sv); printf("Took %lf seconds.\n", (double)(clock() - begin) / CLOCKS_PER_SEC); return 0;
}</lang>
- Output:
The first 50 self numbers are: 1 3 5 7 9 20 31 42 53 64 75 86 97 108 110 121 132 143 154 165 176 187 198 209 211 222 233 244 255 266 277 288 299 310 312 323 334 345 356 367 378 389 400 411 413 424 435 446 457 468 The 100 millionth self number is 1022727208 Took 1.521486 seconds.
Extended
<lang c>#include <stdio.h>
- include <stdlib.h>
- include <time.h>
typedef unsigned char bool;
- define TRUE 1
- define FALSE 0
- define MILLION 1000000LL
- define BILLION 1000 * MILLION
- define MAX_COUNT 103LL*10000*10000 + 11*9 + 1
int digitSum[10000];
void init() {
int i = 9999, s, t, a, b, c, d; for (a = 9; a >= 0; --a) { for (b = 9; b >= 0; --b) { s = a + b; for (c = 9; c >= 0; --c) { t = s + c; for (d = 9; d >= 0; --d) { digitSum[i] = t + d; --i; } } } }
}
void sieve(bool *sv) {
int a, b, c; long long s, n = 0; for (a = 0; a < 103; ++a) { for (b = 0; b < 10000; ++b) { s = digitSum[a] + digitSum[b] + n; for (c = 0; c < 10000; ++c) { sv[digitSum[c]+s] = TRUE; ++s; } n += 10000; } }
}
int main() {
long long count = 0, limit = 1; clock_t begin = clock(), end; bool *p, *sv = (bool*) calloc(MAX_COUNT, sizeof(bool)); init(); sieve(sv); printf("Sieving took %lf seconds.\n", (double)(clock() - begin) / CLOCKS_PER_SEC); printf("\nThe first 50 self numbers are:\n"); for (p = sv; p < sv + MAX_COUNT; ++p) { if (!*p) { if (++count <= 50) { printf("%ld ", p-sv); } else { printf("\n\n Index Self number\n"); break; } } } count = 0; for (p = sv; p < sv + MAX_COUNT; ++p) { if (!*p) { if (++count == limit) { printf("%10lld %11ld\n", count, p-sv); limit *= 10; if (limit == 10 * BILLION) break; } } } free(sv); printf("\nOverall took %lf seconds.\n", (double)(clock() - begin) / CLOCKS_PER_SEC); return 0;
}</lang>
- Output:
Sieving took 7.429969 seconds. The first 50 self numbers are: 1 3 5 7 9 20 31 42 53 64 75 86 97 108 110 121 132 143 154 165 176 187 198 209 211 222 233 244 255 266 277 288 299 310 312 323 334 345 356 367 378 389 400 411 413 424 435 446 457 468 Index Self number 1 1 10 64 100 973 1000 10188 10000 102225 100000 1022675 1000000 10227221 10000000 102272662 100000000 1022727208 1000000000 10227272649 Overall took 11.574158 seconds.
C#
via
(third version) Stripped down, as C# limits the size of an array to Int32.MaxValue, so the sieve isn't large enough to hit the 1,000,000,000th value.
<lang csharp>using System; using static System.Console;
class Program {
const int mc = 103 * 1000 * 10000 + 11 * 9 + 1;
static bool[] sv = new bool[mc + 1];
static void sieve() { int[] dS = new int[10000]; for (int a = 9, i = 9999; a >= 0; a--) for (int b = 9; b >= 0; b--) for (int c = 9, s = a + b; c >= 0; c--) for (int d = 9, t = s + c; d >= 0; d--) dS[i--] = t + d; for (int a = 0, n = 0; a < 103; a++) for (int b = 0, d = dS[a]; b < 1000; b++, n += 10000) for (int c = 0, s = d + dS[b] + n; c < 10000; c++) sv[dS[c] + s++] = true; }
static void Main() { DateTime st = DateTime.Now; sieve(); WriteLine("Sieving took {0}s", (DateTime.Now - st).TotalSeconds); WriteLine("\nThe first 50 self numbers are:"); for (int i = 0, count = 0; count <= 50; i++) if (!sv[i]) { count++; if (count <= 50) Write("{0} ", i); else WriteLine("\n\n Index Self number"); } for (int i = 0, limit = 1, count = 0; i < mc; i++) if (!sv[i]) if (++count == limit) { WriteLine("{0,12:n0} {1,13:n0}", count, i); if (limit == 1e9) break; limit *= 10; } WriteLine("\nOverall took {0}s", (DateTime.Now - st). TotalSeconds); }
}</lang>
- Output:
Timing from tio.run
Sieving took 3.4266187s The first 50 self numbers are: 1 3 5 7 9 20 31 42 53 64 75 86 97 108 110 121 132 143 154 165 176 187 198 209 211 222 233 244 255 266 277 288 299 310 312 323 334 345 356 367 378 389 400 411 413 424 435 446 457 468 Index Self number 1 1 10 64 100 973 1,000 10,188 10,000 102,225 100,000 1,022,675 1,000,000 10,227,221 10,000,000 102,272,662 100,000,000 1,022,727,208 Overall took 7.0237244s
F#
<lang fsharp> // Self numbers. Nigel Galloway: October 6th., 2020 let fN g=let rec fG n g=match n/10 with 0->n+g |i->fG i (g+(n%10)) in fG g g let Self=let rec Self n i g=seq{let g=g@([n..i]|>List.map fN) in yield! List.except g [n..i]; yield! Self (n+100) (i+100) (List.filter(fun n->n>i) g)} in Self 0 99 []
Self |> Seq.take 50 |> Seq.iter(printf "%d "); printfn "" printfn "\n%d" (Seq.item 99999999 Self) </lang>
- Output:
1 3 5 7 9 20 31 42 53 64 75 86 97 108 110 121 132 143 154 165 176 187 198 209 211 222 233 244 255 266 277 288 299 310 312 323 334 345 356 367 378 389 400 411 413 424 435 446 457 468 1022727208
Go
Low memory
Simple-minded, uses very little memory (no sieve) but slow - over 2.5 minutes. <lang go>package main
import (
"fmt" "time"
)
func sumDigits(n int) int {
sum := 0 for n > 0 { sum += n % 10 n /= 10 } return sum
}
func max(x, y int) int {
if x > y { return x } return y
}
func main() {
st := time.Now() count := 0 var selfs []int i := 1 pow := 10 digits := 1 offset := 9 lastSelf := 0 for count < 1e8 { isSelf := true start := max(i-offset, 0) sum := sumDigits(start) for j := start; j < i; j++ { if j+sum == i { isSelf = false break } if (j+1)%10 != 0 { sum++ } else { sum = sumDigits(j + 1) } } if isSelf { count++ lastSelf = i if count <= 50 { selfs = append(selfs, i) if count == 50 { fmt.Println("The first 50 self numbers are:") fmt.Println(selfs) } } } i++ if i%pow == 0 { pow *= 10 digits++ offset = digits * 9 } } fmt.Println("\nThe 100 millionth self number is", lastSelf) fmt.Println("Took", time.Since(st))
}</lang>
- Output:
The first 50 self numbers are: [1 3 5 7 9 20 31 42 53 64 75 86 97 108 110 121 132 143 154 165 176 187 198 209 211 222 233 244 255 266 277 288 299 310 312 323 334 345 356 367 378 389 400 411 413 424 435 446 457 468] The 100 millionth self number is 1022727208 Took 2m35.531949399s
Sieve based
Simple sieve, requires a lot of memory but quick - around 2 seconds.
Nested 'for's used rather than a recursive function for extra speed.
Have also incorporated Enter your username's suggestion (see Talk page) of using partial sums for each loop which improves performance by about 25%. <lang go>package main
import (
"fmt" "time"
)
func sieve() []bool {
sv := make([]bool, 2*1e9+9*9 + 1) n := 0 var s [8]int for a := 0; a < 2; a++ { for b := 0; b < 10; b++ { s[0] = a + b for c := 0; c < 10; c++ { s[1] = s[0] + c for d := 0; d < 10; d++ { s[2] = s[1] + d for e := 0; e < 10; e++ { s[3] = s[2] + e for f := 0; f < 10; f++ { s[4] = s[3] + f for g := 0; g < 10; g++ { s[5] = s[4] + g for h := 0; h < 10; h++ { s[6] = s[5] + h for i := 0; i < 10; i++ { s[7] = s[6] + i for j := 0; j < 10; j++ { sv[s[7]+j+n] = true n++ } } } } } } } } } } return sv
}
func main() {
st := time.Now() sv := sieve() count := 0 fmt.Println("The first 50 self numbers are:") for i := 0; i < len(sv); i++ { if !sv[i] { count++ if count <= 50 { fmt.Printf("%d ", i) } if count == 1e8 { fmt.Println("\n\nThe 100 millionth self number is", i) break } } } fmt.Println("Took", time.Since(st))
}</lang>
- Output:
The first 50 self numbers are: 1 3 5 7 9 20 31 42 53 64 75 86 97 108 110 121 132 143 154 165 176 187 198 209 211 222 233 244 255 266 277 288 299 310 312 323 334 345 356 367 378 389 400 411 413 424 435 446 457 468 The 100 millionth self number is 1022727208 Took 1.984969034s
Extended
This uses horst.h's ideas (see Talk page) to find up to the 1 billionth self number in a reasonable time and using less memory than the simple 'sieve based' approach above would have needed. <lang go>package main
import (
"fmt" "time"
)
const MAX_COUNT = 103*1e4*1e4 + 11*9 + 1
var sv = make([]bool, MAX_COUNT+1) var digitSum = make([]int, 1e4)
func init() {
i := 9999 var s, t int for a := 9; a >= 0; a-- { for b := 9; b >= 0; b-- { s = a + b for c := 9; c >= 0; c-- { t = s + c for d := 9; d >= 0; d-- { digitSum[i] = t + d i-- } } } }
}
func sieve() {
n := 0 for a := 0; a < 103; a++ { for b := 0; b < 1e4; b++ { s := digitSum[a] + digitSum[b] + n for c := 0; c < 1e4; c++ { sv[digitSum[c]+s] = true s++ } n += 1e4 } }
}
func main() {
st := time.Now() sieve() fmt.Println("Sieving took", time.Since(st)) count := 0 fmt.Println("\nThe first 50 self numbers are:") for i := 0; i < len(sv); i++ { if !sv[i] { count++ if count <= 50 { fmt.Printf("%d ", i) } else { fmt.Println("\n\n Index Self number") break } } } count = 0 limit := 1 for i := 0; i < len(sv); i++ { if !sv[i] { count++ if count == limit { fmt.Printf("%10d %11d\n", count, i) limit *= 10 if limit == 1e10 { break } } } } fmt.Println("\nOverall took", time.Since(st))
}</lang>
- Output:
Sieving took 8.286841692s The first 50 self numbers are: 1 3 5 7 9 20 31 42 53 64 75 86 97 108 110 121 132 143 154 165 176 187 198 209 211 222 233 244 255 266 277 288 299 310 312 323 334 345 356 367 378 389 400 411 413 424 435 446 457 468 Index Self number 1 1 10 64 100 973 1000 10188 10000 102225 100000 1022675 1000000 10227221 10000000 102272662 100000000 1022727208 1000000000 10227272649 Overall took 14.647314803s
Julia
The code first bootstraps a sliding window of size 81 and then uses this as a sieve. Note that 81 is the window size because the sum of digits of 999,999,999 (the largest digit sum of a counting number less than 1022727208) is 81. <lang julia>gsum(i) = sum(digits(i)) + i isnonself(i) = any(x -> gsum(x) == i, i-1:-1:i-max(1, ndigits(i)*9)) const last81 = filter(isnonself, 1:5000)[1:81]
function checkselfnumbers()
i, selfcount = 1, 0 while selfcount <= 100_000_000 && i <= 1022727208 if !(i in last81) selfcount += 1 if selfcount < 51 print(i, " ") elseif selfcount == 51 println() elseif selfcount == 100_000_000 println(i == 1022727208 ? "Yes, $i is the 100,000,000th self number." : "No, instead $i is the 100,000,000th self number.") end end popfirst!(last81) push!(last81, gsum(i)) i += 1 end
end
checkselfnumbers()
</lang>
- Output:
1 3 5 7 9 20 31 42 53 64 75 86 97 108 110 121 132 143 154 165 176 187 198 209 211 222 233 244 255 266 277 288 299 310 312 323 334 345 356 367 378 389 400 411 413 424 435 446 457 468 Yes, 1022727208 is the 100,000,000th self number.
Faster version
Contains tweaks peculiar to the "10 to the nth" self number. Timings include compilation times. <lang julia>const MAXCOUNT = 103 * 10000 * 10000 + 11 * 9 + 1
function dosieve!(sieve, digitsum9999)
n = 1 for a in 1:103, b in 1:10000 s = digitsum9999[a] + digitsum9999[b] + n for c in 1:10000 sieve[digitsum9999[c] + s] = true s += 1 end n += 10000 end
end
initdigitsum() = reverse!(vec([sum(k) for k in Iterators.product(9:-1:0, 9:-1:0, 9:-1:0, 9:-1:0)]))
function findselves()
sieve = zeros(Bool, MAXCOUNT+1) println("Sieve time:") @time begin digitsum = initdigitsum() dosieve!(sieve, digitsum) end cnt = 1 for i in 1:MAXCOUNT+1 if !sieve[i] cnt > 50 && break print(i, " ") cnt += 1 end end println() limit, cnt = 1, 0 for i in 0:MAXCOUNT cnt += 1 - sieve[i + 1] if cnt == limit println(lpad(cnt, 10), lpad(i, 12)) limit *= 10 end end
end
@time findselves()
</lang>
- Output:
Sieve time: 7.187635 seconds (2 allocations: 78.203 KiB) 1 3 5 7 9 20 31 42 53 64 75 86 97 108 110 121 132 143 154 165 176 187 198 209 211 222 233 244 255 266 277 288 299 310 312 323 334 345 356 367 378 389 400 411 413 424 435 446 457 468 1 1 10 64 100 973 1000 10188 10000 102225 100000 1022675 1000000 10227221 10000000 102272662 100000000 1022727208 1000000000 10227272649 16.999383 seconds (42.92 k allocations: 9.595 GiB, 0.01% gc time)
Pascal
Just "sieving" with only one follower of every number
Extended to 10.23e9 <lang pascal>program selfnumb; {$IFDEF FPC}
{$MODE Delphi} {$Optimization ON,ALL}
{$IFEND} {$IFDEF DELPHI} {$APPTYPE CONSOLE} {$IFEND} uses
sysutils;
const
MAXCOUNT =103*10000*10000+11*9+ 1;
type
tDigitSum9999 = array[0..9999] of Uint8; tpDigitSum9999 = ^tDigitSum9999;
var
DigitSum9999 : tDigitSum9999; sieve : array of boolean;
procedure dosieve; var
pSieve : pBoolean; pDigitSum :tpDigitSum9999; n,c,b,a,s : NativeInt;
Begin
pSieve := @sieve[0]; pDigitSum := @DigitSum9999[0]; n := 0; for a := 0 to 102 do for b := 0 to 9999 do Begin s := pDigitSum^[a]+pDigitSum^[b]+n; for c := 0 to 9999 do Begin pSieve[pDigitSum^[c]+s] := true; s+=1; end; inc(n,10000); end;
end;
procedure InitDigitSum; var
i,d,c,b,a : NativeInt;
begin
i := 9999; for a := 9 downto 0 do for b := 9 downto 0 do for c := 9 downto 0 do for d := 9 downto 0 do Begin DigitSum9999[i] := a+b+c+d; dec(i); end;
end;
procedure OutPut(cnt,i:NativeUint); Begin
writeln(cnt:10,i:12);
end;
var
pSieve : pboolean; T0 : Uint64; i,cnt,limit,One: NativeUInt;
BEGIN
setlength(sieve,MAXCOUNT); pSieve := @sieve[0]; T0 := GetTickCount64; InitDigitSum; dosieve; writeln('Sievetime : ',(GetTickCount64-T0 )/1000:8:3,' sec'); //find first 50 cnt := 0; for i := 0 to MAXCOUNT do Begin if NOT(pSieve[i]) then Begin inc(cnt); if cnt <= 50 then write(i:4) else BREAK; end; end; writeln; One := 1; limit := One; cnt := 0; for i := 0 to MAXCOUNT do Begin inc(cnt,One-Ord(pSieve[i])); if cnt = limit then Begin OutPut(cnt,i); limit := limit*10; end; end;
END.</lang>
- Output:
time ./selfnumb Sievetime : 6.579 sec 1 3 5 7 9 20 31 42 53 64 75 86 97 108 110 121 132 143 154 165 176 187 198 209 211 222 233 244 255 266 277 288 299 310 312 323 334 345 356 367 378 389 400 411 413 424 435 446 457 468 1 1 10 64 100 973 1000 10188 10000 102225 100000 1022675 1000000 10227221 10000000 102272662 100000000 1022727208 1000000000 10227272649 real 0m13,252s
Phix
Replacing the problematic sv[a+b+... line with a bit of dirty inline assembly improved performance by 90%
(Of course you lose bounds checking, type checking, negative subscripts, fraction handling, and all that jazz.)
We use a string of Y/N for the sieve to force one byte per element ('\0' and 1 would be equally valid).
<lang Phix>if machine_bits()=32 then crash("requires 64 bit") end if
function sieve()
string sv = repeat('N',2*1e9+9*9+1) integer n = 0 for a=0 to 1 do for b=0 to 9 do for c=0 to 9 do for d=0 to 9 do for e=0 to 9 do for f=0 to 9 do for g=0 to 9 do for h=0 to 9 do for i=0 to 9 do for j=0 to 9 do
-- n += 1 -- sv[a+b+c+d+e+f+g+h+i+j+n] = 'Y'
- ilASM{
[32] -- (allows clean compilation on 32 bit, before crash as above) [64] mov rax,[a] mov r12,[b] mov r13,[c] mov r14,[d] mov r15,[e] add r12,r13 add r14,r15 add rax,r12 mov rdi,[sv] add rax,r14 mov r12,[f] mov r13,[g] mov r14,[h] mov r15,[i] add r12,r13 shl rdi,2 mov rcx,[n] mov r13,[j] add r14,r15 add rax,r12 add rax,r14 add r13,rcx add rax,r13 add rcx,1 mov byte[rdi+rax],'Y' mov [n],rcx } end for end for end for end for end for end for end for end for printf(1,"%d,%d\r",{a,b}) -- (show progress) end for end for return sv
end function
atom t0 = time() string sv = sieve() printf(1,"sieve build took %s\n",{elapsed(time()-t0)}) integer count = 0 printf(1,"The first 50 self numbers are:\n") for i=1 to length(sv) do
if sv[i]='N' then count += 1 if count <= 50 then printf(1,"%3d ",i-1) if remainder(count,10)=0 then printf(1,"\n") end if end if if count == 1e8 then printf(1,"\nThe %,dth self number is %,d (%s)\n", {count,i-1,elapsed(time()-t0)}) exit end if end if
end for</lang>
- Output:
sieve build took 6.6s The first 50 self numbers are: 1 3 5 7 9 20 31 42 53 64 75 86 97 108 110 121 132 143 154 165 176 187 198 209 211 222 233 244 255 266 277 288 299 310 312 323 334 345 356 367 378 389 400 411 413 424 435 446 457 468 The 100,000,000th self number is 1,022,727,208 (11.0s)
generator dictionary
While this is dog-slow (see shocking estimate below), it is interesting to note that even by the time it generates the
10,000,000th number, it is only having to juggle a mere 27 generators.
Just a shame that we had to push over 10,000,000 generators onto that stack, and tried to push quite a few more.
Memory use is pretty low, around ~4MB.
[unlike the above, this is perfectly happy on both 32 and 64 bit]
Long story short: this works much the same as a prime sieve, in which you only need to eliminate multiples
of previous primes. Here, you only need to eliminate digital additions of previous safe numbers.
Only after writing this did I understand how to write a sliding sieve, which it turns out is a much better idea (see below).
Still, this is pretty interesting and quite neat.
<lang Phix>integer gd = new_dict(), gdhead = 2, n = 0
function ng(integer n)
integer r = n while r do n += remainder(r,10) r = floor(r/10) end while return n
end function
function self(integer /*i*/) -- note: assumes sequentially invoked (arg i unused)
n += 1 while n=gdhead do gdhead = pop_dict(gd)[1] setd(ng(gdhead),0,gd) n += (n!=gdhead) end while setd(ng(n),0,gd) return n
end function
atom t0 = time() printf(1,"The first 50 self numbers are:\n") pp(apply(tagset(50),self),{pp_IntFmt,"%3d",pp_IntCh,false})
constant limit = 10000000 integer chk = 100 printf(1,"\n i n size time\n") for i=51 to limit do
n = self(i) if i=chk then printf(1,"%,11d %,11d %6d %s\n",{i,n,dict_size(gd),elapsed(time()-t0)}) chk *= 10 end if
end for printf(1,"\nEstimated time for %,d :%s\n",{1e8,elapsed((time()-t0)*1e8/limit)})</lang>
- Output:
The first 50 self numbers are: { 1, 3, 5, 7, 9, 20, 31, 42, 53, 64, 75, 86, 97,108,110,121,132,143, 154,165,176,187,198,209,211,222,233,244,255,266,277,288,299,310,312,323, 334,345,356,367,378,389,400,411,413,424,435,446,457,468} i n size time 100 973 18 0.1s 1,000 10,188 13 0.2s 10,000 102,225 10 1.0s 100,000 1,022,675 20 9.3s 1,000,000 10,227,221 17 1 minute and 37s 10,000,000 102,272,662 27 16 minutes and 40s Estimated time for 100,000,000 :2 hours, 46 minutes and 37s
For the record, I would not be at all surprised should a translation of this beat 20 minutes (for 1e8)
sliding sieve
Mid-speed, perhaps helped by a slightly smarter way of calculating/updating the digit sums
Similar to some other entries, we (only) need a sieve of 9*9, +1 here as I test an entry after the slide.
<lang Phix>
--sequence sieve = repeat(0,82), -- (~25% slower)
sequence sieve = repeat(0,8192),
digits = repeat(0,10)
integer offset = 0,
digit_sum = 0, n = 0
procedure next_self()
while true do n += 1 for i=length(digits) to 1 by -1 do integer d = digits[i] if d!=9 then digits[i] = d+1 digit_sum += 1 exit end if digits[i] = 0 digit_sum -= 9 end for integer k = n+digit_sum-offset if k>length(sieve) then integer j = 1 for i=n-offset to length(sieve) do sieve[j] = sieve[i] j += 1 end for sieve[j..$] = 0 offset = n-1 k = digit_sum+1 end if sieve[k] = 1 if sieve[n-offset]=0 then exit end if end while -- (result in n)
end procedure
constant limit = 100000000 atom t0 = time() printf(1,"The first 50 self numbers are:\n") integer chk = 100 for i=1 to limit do
next_self() if i<=50 then printf(1," %3d%s",{n,iff(mod(i,25)=0?"\n":"")}) elsif i=chk then if chk=100 then printf(1,"\n i n time\n") end if printf(1,"%,12d %,13d %s\n",{i,n,elapsed(time()-t0)}) chk *= 10 end if
end for printf(1,"\nEstimated time for %,d :%s\n",{1e9,elapsed((time()-t0)*1e9/limit)})</lang>
- Output:
The first 50 self numbers are: 1 3 5 7 9 20 31 42 53 64 75 86 97 108 110 121 132 143 154 165 176 187 198 209 211 222 233 244 255 266 277 288 299 310 312 323 334 345 356 367 378 389 400 411 413 424 435 446 457 468 i n time 100 973 0.1s 1,000 10,188 0.1s 10,000 102,225 0.1s 100,000 1,022,675 0.1s 1,000,000 10,227,221 0.5s 10,000,000 102,272,662 4.5s 100,000,000 1,022,727,208 44.5s Estimated time for 1,000,000,000 :7 minutes and 25s
REXX
first 50 self numbers
<lang rexx>/*REXX program displays N self numbers, aka Colombian or Devlali numbers. OEIS A3052.*/ parse arg n . /*obtain optional argument from the CL.*/ if n== | n=="," then n= 50 /*Not specified? Then use the default.*/ @.=. /*initialize the array of self numbers.*/
do j=1 for n*10 /*scan through ten times the #s wanted.*/ $= j /*1st part of sum is the number itself.*/ do k=1 for length(j) /*sum the decimal digits in the number.*/ $= $ + substr(j, k, 1) /*add a particular digit to the sum. */ end /*k*/ @.$= /*mark J as not being a self number. */ end /*j*/
list= 1 /*initialize list to the first number. */
- = 1 /*the count of self numbers (so far). */
do i=3 until #==n /*traipse through array 'til satisfied.*/ if @.i== then iterate /*Not a self number? Then skip it. */ #= # + 1 /*bump the counter of self numbers. */ list= list i /*append it to the list of self numbers*/ end /*i*/
say 'the first ' n " self numbers are:" /*display the title for the output list*/ say list /*display list of self numbers ──►term.*/ exit 0 /*stick a fork in it, we're all done. */</lang>
- output when using the default input:
the first 50 self numbers are: 1 3 5 7 9 20 31 42 53 64 75 86 97 108 110 121 132 143 154 165 176 187 198 209 211 222 233 244 255 266 277 288 299 310 312 323 334 345 356 367 378 389 400 411 413 424 435 446 457 468
ten millionth self number
<lang rexx>/*REXX pgm displays the Nth self number, aka Colombian or Devlali numbers. OEIS A3052.*/ numeric digits 20 /*ensure enough decimal digits for #. */ parse arg high . /*obtain optional argument from the CL.*/ if high== | high=="," then high= 100000000 /*Not specified? Then use 100M default*/ i= 1; pow= 10; digs= 1; offset= 9; $= 0 /*$: the last self number found. */
- = 0 /*count of self numbers (so far). */
do while #<high; isSelf= 1 /*assume a self number (so far). */ start= max(i-offset, 0) sum= sumDigs(start)
do j=start to i-1 if j+sum==i then do; isSelf= 0 /*found a non self number. */ iterate /*keep looking for more self numbers. */ end jp= j + 1 /*shortcut variable for next statement.*/ if jp//10==0 then sum= sumDigs(jp) else sum= sum + 1 end /*j*/
if isSelf then do; #= # + 1 /*bump the count of self numbers. */ $= i /*save the last self number found. */ end i= i + 1 if i//pow==0 then do; pow= pow * 10 digs= digs + 1 /*bump the number of decimal digits. */ offset= digs * 9 /*bump the offset by a factor of nine. */ end end /*while*/
say say 'the ' commas(high)th(high) " self number is: " commas($) exit 0 /*stick a fork in it, we're all done. */ /*──────────────────────────────────────────────────────────────────────────────────────*/ sumDigs: parse arg s 2 x; do k=1 for length(x) /*get 1st dig, & also get the rest.*/
s= s + substr(x, k, 1) /*add a particular digit to the sum.*/ end /*k*/; return s
/*──────────────────────────────────────────────────────────────────────────────────────*/ commas: parse arg _; do c=length(_)-3 to 1 by -3; _=insert(',', _, c); end; return _ th: parse arg th; return word('th st nd rd',1+(th//10)*(th//100%10\==1)*(th//10<4))</lang>
- output when using the default input:
the 100,000,000th self number is: 1,022,727,208
Wren
Just the sieve based version as the low memory version would take too long to run in Wren.
Note that you need a lot of memory to run this as Bools in Wren require 8 bytes of storage compared to 1 byte in Go.
Unsurprisingly, very slow compared to the Go version as Wren is interpreted and uses floating point arithmetic for all numerical work. <lang ecmascript>var sieve = Fn.new {
var sv = List.filled(2*1e9+9*9+1, false) var n = 0 var s = [0] * 8 for (a in 0..1) { for (b in 0..9) { s[0] = a + b for (c in 0..9) { s[1] = s[0] + c for (d in 0..9) { s[2] = s[1] + d for (e in 0..9) { s[3] = s[2] + e for (f in 0..9) { s[4] = s[3] + f for (g in 0..9) { s[5] = s[4] + g for (h in 0..9) { s[6] = s[5] + h for (i in 0..9) { s[7] = s[6] + i for (j in 0..9) { sv[s[7] + j + n] = true n = n + 1 } } } } } } } } } } return sv
}
var st = System.clock var sv = sieve.call() var count = 0 System.print("The first 50 self numbers are:") for (i in 0...sv.count) {
if (!sv[i]) { count = count + 1 if (count <= 50) System.write("%(i) ") if (count == 1e8) { System.print("\n\nThe 100 millionth self number is %(i)") break } }
} System.print("Took %(System.clock-st) seconds.")</lang>
- Output:
The first 50 self numbers are: 1 3 5 7 9 20 31 42 53 64 75 86 97 108 110 121 132 143 154 165 176 187 198 209 211 222 233 244 255 266 277 288 299 310 312 323 334 345 356 367 378 389 400 411 413 424 435 446 457 468 The 100 millionth self number is 1022727208 Took 222.789713 seconds.