Self numbers: Difference between revisions

From Rosetta Code
Content added Content deleted
m (→‎sliding sieve: removed blank line)
(→‎{{header|AppleScript}}: Several modifications for simplicity, efficiency, and hopefully clarity. Added a "fast forward" function which reduces the time taken to get the 100,000,000th number from just over two minutes to less than a millisecond.)
Line 15: Line 15:
=={{header|AppleScript}}==
=={{header|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.
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.


<lang applescript>(*
<lang applescript>(*
Base-10 self numbers by index (single or range).
Base-10 self numbers by index (single or range).
After the initial single-digit numbers, the sequence follows an observed pattern whereby self numbers
Follows an observed sequence pattern whereby, after the initial single-digit odd numbers, self numbers are
occur in groups whose consecutive members are 11 apart. The interval between these groups is 2.
grouped in runs whose members occur at numeric intervals of 11. Runs after the first one come in blocks of
They occur in blocks of ten: eight ten-number groups sandwiched between two shorter groups.
ten: eight runs of ten numbers followed by two shorter runs. The numeric interval between runs is usually 2,
The intervals between blocks and the lengths of the shorter groups depend the most signficant
but that between shorter runs, and their length, depend on the highest-order digit change occurring in them.
This connection with significant digit change means every ten blocks form a higher-order block, every ten
digit to have changed in the number sequence since they were last set.
of these a higher-order-still block, and so on.
*)
*)
on selfNumbers(indexRange)
on selfNumbers(indexRange)
-- Subscript to monitor progress, store required results, and indicate when finished.
set indexRange to indexRange as list
set indexRange to indexRange as list
-- Script object with subhandlers and associated properties.
script storer
script |subscript|
property startIndex : beginning of indexRange
property startIndex : beginning of indexRange
property endIndex : end of indexRange
property endIndex : end of indexRange
property counter : 0
property counter : 0
property currentSelf : -1
property output : {}
property output : {}
-- Advance to the next self number in the sequence, append it to the output if required, indicate if finished.
on doneAfter(thisSelf)
on doneAfterAdding(interval)
set currentSelf to currentSelf + interval
set counter to counter + 1
set counter to counter + 1
if (counter < startIndex) then return false
if (counter < startIndex) then return false
set end of my output to thisSelf
set end of my output to currentSelf
return (counter = endIndex)
return (counter = endIndex)
end doneAfter
end doneAfterAdding
-- If necessary, fast forward to the last self number before the lowest-order block containing the first number required.
on fastForward()
if (counter startIndex) then return
-- The highest-order blocks whose ends this script is thought to handle correctly contain 9,777,777,778 self numbers.
-- The difference between equivalently positioned numbers in these blocks is 100,000,000,001.
-- The figures for successively lower-order blocks are obtained by successively removing 7s and 0s!
set indexDiff to 9.777777778E+9
set numericDiff to 1.00000000001E+11
repeat until ((indexDiff < 98) or (counter = startIndex))
set test to counter + indexDiff
if (test > startIndex) then
set indexDiff to indexDiff div 10 + 1
set numericDiff to numericDiff div 10 + 1
else
set counter to test
set currentSelf to (currentSelf + numericDiff)
end if
end repeat
end fastForward
-- Work out a shorter run length based on the most significant digit change about to happen.
on getShorterRunLength()
set shorterRunLength to 9
set temp to (|subscript|'s currentSelf) div 1000
repeat while (temp mod 10 is 9)
set shorterRunLength to shorterRunLength - 1
set temp to temp div 10
end repeat
return shorterRunLength
end getShorterRunLength
end script
end script
-- Start with the single-digit odds.
-- Main process. Start with the single-digit odd numbers and first run.
set thisSelf to -1
repeat 5 times
if (|subscript|'s doneAfterAdding(2)) then return |subscript|'s output
repeat 4 times
set thisSelf to thisSelf + 2
if (storer's doneAfter(thisSelf)) then return storer's output
end repeat
end repeat
repeat 9 times
if (|subscript|'s doneAfterAdding(11)) then return |subscript|'s output
end repeat
-- Fast forward if the start index hasn't yet been reached.
tell |subscript| to fastForward()
-- Main loop.
-- Sequencing loop, per lowest-order block.
repeat
set shortGroupLength to 10
-- Eight ten-number runs, each at a numeric interval of 2 from the end of the previous one.
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
repeat 8 times
if (|subscript|'s doneAfterAdding(2)) then return |subscript|'s output
set thisSelf to thisSelf + 2
if (storer's doneAfter(thisSelf)) then return storer's output
repeat 9 times
repeat 9 times
set thisSelf to thisSelf + 11
if (|subscript|'s doneAfterAdding(11)) then return |subscript|'s output
if (storer's doneAfter(thisSelf)) then return storer's output
end repeat
end repeat
end repeat
end repeat
-- Two shorter runs, the second at an interval inversely related to their length.
set shorterRunLength to |subscript|'s getShorterRunLength()
-- Second short group. Its length depends on the most significant digit to tick over
repeat with interval in {2, 2 + (10 - shorterRunLength) * 13}
-- at the thousands boundary within it. The interval preceding it is 2.
if (|subscript|'s doneAfterAdding(interval)) then return |subscript|'s output
set thisSelf to thisSelf + 2
repeat (shorterRunLength - 1) times
if (storer's doneAfter(thisSelf)) then return storer's output
if (|subscript|'s doneAfterAdding(11)) then return |subscript|'s output
set i to 2
set shortGroupLength to 9
end repeat
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
end repeat
end repeat
return storer's output
return |subscript|'s output
end selfNumbers
end selfNumbers


Line 102: Line 116:
-- One hundred millionth:
-- One hundred millionth:
selfNumbers(100000000)
selfNumbers(100000000)

--> {1022727208}</lang>
--> {1.022727208E+9}</lang>


=={{header|C}}==
=={{header|C}}==

Revision as of 11:04, 19 October 2020

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

<lang applescript>(*

   Base-10 self numbers by index (single or range).
   Follows an observed sequence pattern whereby, after the initial single-digit odd numbers, self numbers are
   grouped in runs whose members occur at numeric intervals of 11. Runs after the first one come in blocks of
   ten: eight runs of ten numbers followed by two shorter runs. The numeric interval between runs is usually 2,
   but that between shorter runs, and their length, depend on the highest-order digit change occurring in them.
   This connection with significant digit change means every ten blocks form a higher-order block, every ten
   of these a higher-order-still block, and so on.
  • )

on selfNumbers(indexRange)

   set indexRange to indexRange as list
   -- Script object with subhandlers and associated properties.
   script |subscript|
       property startIndex : beginning of indexRange
       property endIndex : end of indexRange
       property counter : 0
       property currentSelf : -1
       property output : {}
       
       -- Advance to the next self number in the sequence, append it to the output if required, indicate if finished.
       on doneAfterAdding(interval)
           set currentSelf to currentSelf + interval
           set counter to counter + 1
           if (counter < startIndex) then return false
           set end of my output to currentSelf
           return (counter = endIndex)
       end doneAfterAdding
       
       -- If necessary, fast forward to the last self number before the lowest-order block containing the first number required.
       on fastForward()
           if (counter ≥ startIndex) then return
           -- The highest-order blocks whose ends this script is thought to handle correctly contain 9,777,777,778 self numbers.
           -- The difference between equivalently positioned numbers in these blocks is 100,000,000,001.
           -- The figures for successively lower-order blocks are obtained by successively removing 7s and 0s!
           set indexDiff to 9.777777778E+9
           set numericDiff to 1.00000000001E+11
           repeat until ((indexDiff < 98) or (counter = startIndex))
               set test to counter + indexDiff
               if (test > startIndex) then
                   set indexDiff to indexDiff div 10 + 1
                   set numericDiff to numericDiff div 10 + 1
               else
                   set counter to test
                   set currentSelf to (currentSelf + numericDiff)
               end if
           end repeat
       end fastForward
       
       -- Work out a shorter run length based on the most significant digit change about to happen.
       on getShorterRunLength()
           set shorterRunLength to 9
           set temp to (|subscript|'s currentSelf) div 1000
           repeat while (temp mod 10 is 9)
               set shorterRunLength to shorterRunLength - 1
               set temp to temp div 10
           end repeat
           return shorterRunLength
       end getShorterRunLength
   end script
   
   -- Main process. Start with the single-digit odd numbers and first run.
   repeat 5 times
       if (|subscript|'s doneAfterAdding(2)) then return |subscript|'s output
   end repeat
   repeat 9 times
       if (|subscript|'s doneAfterAdding(11)) then return |subscript|'s output
   end repeat
   -- Fast forward if the start index hasn't yet been reached.
   tell |subscript| to fastForward()
   
   -- Sequencing loop, per lowest-order block.
   repeat
       -- Eight ten-number runs, each at a numeric interval of 2 from the end of the previous one.
       repeat 8 times
           if (|subscript|'s doneAfterAdding(2)) then return |subscript|'s output
           repeat 9 times
               if (|subscript|'s doneAfterAdding(11)) then return |subscript|'s output
           end repeat
       end repeat
       -- Two shorter runs, the second at an interval inversely related to their length.
       set shorterRunLength to |subscript|'s getShorterRunLength()
       repeat with interval in {2, 2 + (10 - shorterRunLength) * 13}
           if (|subscript|'s doneAfterAdding(interval)) then return |subscript|'s output
           repeat (shorterRunLength - 1) times
               if (|subscript|'s doneAfterAdding(11)) then return |subscript|'s output
           end repeat
       end repeat
   end repeat
   
   return |subscript|'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)

--> {1.022727208E+9}</lang>

C

Sieve based

Translation of: Go

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>

  1. include <stdlib.h>
  2. include <time.h>

typedef unsigned char bool;

  1. define TRUE 1
  2. define FALSE 0
  3. define MILLION 1000000
  4. define BILLION 1000 * MILLION
  5. 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

Translation of: Pascal

<lang c>#include <stdio.h>

  1. include <stdlib.h>
  2. include <time.h>

typedef unsigned char bool;

  1. define TRUE 1
  2. define FALSE 0
  3. define MILLION 1000000LL
  4. define BILLION 1000 * MILLION
  5. 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#

Translation of: Pascal

via

Translation of: Go

(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

Translation of: Pascal

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

Translation of: Pascal

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

Works with: Free Pascal


Just "sieving" with only one follower of every number

Translation of: Go

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

Translation of: Go

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) -- (1.86GB)
   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'

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

Translation of: Go (low memory)

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

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

Translation of: Go

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.