I'm working on modernizing Rosetta Code's infrastructure. Starting with communications. Please accept this time-limited open invite to RC's Slack.. --Michael Mol (talk) 20:59, 30 May 2020 (UTC)

Untouchable numbers

From Rosetta Code
Task
Untouchable numbers
You are encouraged to solve this task according to the task description, using any language you may know.
Definitions
  •   Untouchable numbers   are also known as   nonaliquot numbers.
  •   An   untouchable number   is a positive integer that cannot be expressed as the sum of all the proper divisors of any positive integer.   (From Wikipedia)
  •   The   sum of all the proper divisors   is also known as   the   aliquot sum.
  •   An   untouchable   are those numbers that are not in the image of the aliquot sum function.   (From Wikipedia)
  •   Untouchable numbers:   impossible values for the sum of all aliquot parts function.   (From OEIS:   The On-line Encyclopedia of Integer Sequences®)
  •   An untouchable number is a positive integer that is not the sum of the proper divisors of any number.   (From MathWorld™)


Observations and conjectures

All untouchable numbers > 5 are composite numbers.

No untouchable number is perfect.

No untouchable number is sociable.

No untouchable number is a Mersenne prime.

No untouchable number is   one more   than a prime number,   since if   p   is prime,   then the sum of the proper divisors of   p2   is  p + 1.

No untouchable number is   three more   than an odd prime number,   since if   p   is an odd prime,   then the sum of the proper divisors of   2p   is  p + 3.

The number  5  is believed to be the only odd untouchable number,   but this has not been proven:   it would follow from a slightly stronger version of the   Goldbach's conjecture,   since the sum of the proper divisors of   pq   (with   p, q   being distinct primes)   is   1 + p + q.

There are infinitely many untouchable numbers,   a fact that was proven by   Paul Erdős.

According to Chen & Zhao,   their natural density is at least   d > 0.06.


Task
  •   show  (in a grid format)  all untouchable numbers  ≤  2,000.
  •   show (for the above)   the   count   of untouchable numbers.
  •   show the   count   of untouchable numbers from unity up to   (inclusive):
  •                   10
  •                 100
  •               1,000
  •             10,000
  •           100,000
  •   ... or as high as is you think is practical.
  •   all output is to be shown here, on this page.


See also



C++[edit]

This solution implements Talk:Untouchable_numbers#Nice_recursive_solution

The Function[edit]

 
// Untouchable Numbers : Nigel Galloway - March 4th., 2021;
#include <functional>
#include <bitset>
#include <iostream>
#include <cmath>
using namespace std; using Z0=long long; using Z1=optional<Z0>; using Z2=optional<array<int,3>>; using Z3=function<Z2()>;
const int maxUT{3000000}, dL{(int)log2(maxUT)};
struct uT{
bitset<maxUT+1>N; vector<int> G{}; array<Z3,int(dL+1)>L{Z3{}}; int sG{0},mUT{};
void _g(int n,int g){if(g<=mUT){N[g]=false; return _g(n,n+g);}}
Z1 nxt(const int n){if(n>mUT) return Z1{}; if(N[n]) return Z1(n); return nxt(n+1);}
Z3 fN(const Z0 n,const Z0 i,int g){return [=]()mutable{if(g<sG && ((n+i)*(1+G[g])-n*G[g]<=mUT)) return Z2{{n,i,g++}}; return Z2{};};}
Z3 fG(Z0 n,Z0 i,const int g){Z0 e{n+i},l{1},p{1}; return [=]()mutable{n=n*G[g]; p=p*G[g]; l=l+p; i=e*l-n; if(i<=mUT) return Z2{{n,i,g}}; return Z2{};};}
void fL(Z3 n, int g){for(;;){
if(auto i=n()){N[(*i)[1]]=false; L[g+1]=fN((*i)[0],(*i)[1],(*i)[2]+1); g=g+1; continue;}
if(auto i=L[g]()){n=fG((*i)[0],(*i)[1],(*i)[2]); continue;}
if(g>0) if(auto i=L[g-1]()){ g=g-1; n=fG((*i)[0],(*i)[1],(*i)[2]); continue;}
if(g>0){ n=[](){return Z2{};}; g=g-1; continue;} break;}
}
int count(){int g{0}; for(auto n=nxt(0); n; n=nxt(*n+1)) ++g; return g;}
uT(const int n):mUT{n}{
N.set(); N[0]=false; N[1]=false; for(auto n=nxt(0);*n<=sqrt(mUT);n=nxt(*n+1)) _g(*n,*n+*n); for(auto n=nxt(0); n; n=nxt(*n+1)) G.push_back(*n); sG=G.size();
N.set(); N[0]=false; L[0]=fN(1,0,0); fL([](){return Z2{};},0);
}
};
 

The Task[edit]

Less than 2000
 
int main(int argc, char *argv[]) {
int c{0}; auto n{uT{2000}}; for(auto g=n.nxt(0); g; g=n.nxt(*g+1)){if(c++==30){c=1; printf("\n");} printf("%4d ",*g);} printf("\n");
}
 
Output:
   2    5   52   88   96  120  124  146  162  188  206  210  216  238  246  248  262  268  276  288  290  292  304  306  322  324  326  336  342  372
 406  408  426  430  448  472  474  498  516  518  520  530  540  552  556  562  576  584  612  624  626  628  658  668  670  708  714  718  726  732
 738  748  750  756  766  768  782  784  792  802  804  818  836  848  852  872  892  894  896  898  902  926  934  936  964  966  976  982  996 1002
1028 1044 1046 1060 1068 1074 1078 1080 1102 1116 1128 1134 1146 1148 1150 1160 1162 1168 1180 1186 1192 1200 1212 1222 1236 1246 1248 1254 1256 1258
1266 1272 1288 1296 1312 1314 1316 1318 1326 1332 1342 1346 1348 1360 1380 1388 1398 1404 1406 1418 1420 1422 1438 1476 1506 1508 1510 1522 1528 1538
1542 1566 1578 1588 1596 1632 1642 1650 1680 1682 1692 1716 1718 1728 1732 1746 1758 1766 1774 1776 1806 1816 1820 1822 1830 1838 1840 1842 1844 1852
1860 1866 1884 1888 1894 1896 1920 1922 1944 1956 1958 1960 1962 1972 1986 1992
Count less than or equal 100000
 
int main(int argc, char *argv[]) {
int z{100000}; auto n{uT{z}}; cout<<"untouchables below "<<z<<"->"<<n.count()<<endl;
}
 
Output:
untouchables below 100000->13863
real    0m2.928s
Count less than or equal 1000000
 
int main(int argc, char *argv[]) {
int z{1000000}; auto n{uT{z}}; cout<<"untouchables below "<<z<<"->"<<n.count()<<endl;
}
 
Output:
untouchables below 1000000->150232
real    3m6.909s
Count less than or equal 2000000
 
int main(int argc, char *argv[]) {
int z{2000000}; auto n{uT{z}}; cout<<"untouchables below "<<z<<"->"<<n.count()<<endl;
}
 
Output:
untouchables below 2000000->305290
real    11m28.031s

Delphi[edit]

Translation of: Go
 
program Untouchable_numbers;
 
{$APPTYPE CONSOLE}
 
uses
System.SysUtils;
 
function SumDivisors(n: Integer): Integer;
begin
Result := 1;
var k := 2;
if not odd(n) then
k := 1;
var i := 1 + k;
while i * i <= n do
begin
if (n mod i) = 0 then
begin
inc(Result, i);
var j := n div i;
if j <> i then
inc(Result, j);
end;
inc(i, k);
end;
end;
 
function Sieve(n: Integer): TArray<Boolean>;
begin
inc(n);
SetLength(result, n + 1);
for var i := 6 to n do
begin
var sd := SumDivisors(i);
if sd <= n then
result[sd] := True;
end;
end;
 
function PrimeSieve(limit: Integer): TArray<Boolean>;
begin
inc(limit);
SetLength(result, limit);
Result[0] := True;
Result[1] := True;
 
var p := 3;
repeat
var p2 := p * p;
if p2 >= limit then
Break;
var i := p2;
while i < limit do
begin
 
Result[i] := True;
inc(i, 2 * p);
end;
 
repeat
inc(p, 2);
until not Result[p];
 
until (False);
 
end;
 
function Commatize(n: Double): string;
var
fmt: TFormatSettings;
begin
fmt := TFormatSettings.Create('en-US');
Result := n.ToString(ffNumber, 64, 0, fmt);
end;
 
begin
var limit := 1000000;
var c := primeSieve(limit);
var s := sieve(63 * limit);
var untouchable: TArray<Integer> := [2, 5];
var n := 6;
while n <= limit do
begin
if not s[n] and c[n - 1] and c[n - 3] then
begin
SetLength(untouchable, Length(untouchable) + 1);
untouchable[High(untouchable)] := n;
end;
inc(n, 2);
end;
writeln('List of untouchable numbers <= 2,000:');
var count := 0;
var i := 0;
while untouchable[i] <= 2000 do
begin
write(commatize(untouchable[i]): 6);
if ((i + 1) mod 10) = 0 then
writeln;
inc(i);
end;
writeln(#10#10, commatize(count): 7, ' untouchable numbers were found <= 2,000');
 
var p := 10;
count := 0;
for n in untouchable do
begin
inc(count);
if n > p then
begin
var cc := commatize(count - 1);
var cp := commatize(p);
writeln(cc, ' untouchable numbers were found <= ', cp);
p := p * 10;
if p = limit then
Break;
end;
end;
 
var cu := commatize(Length(untouchable));
var cl := commatize(limit);
writeln(cu:7, ' untouchable numbers were found <= ', cl);
readln;
end.

F#[edit]

The Function[edit]

This task uses Extensible Prime Generator (F#).
It implements Talk:Untouchable_numbers#Nice_recursive_solution

 
// Applied dendrology. Nigel Galloway: February 15., 2021
let uT a=let N,G=Array.create(a+1) true, [|yield! primes64()|>Seq.takeWhile((>)(int64 a))|]
let fN n i e=let mutable p=e-1 in (fun()->p<-p+1; if p<G.Length && (n+i)*(1L+G.[p])-n*G.[p]<=(int64 a) then Some(n,i,p) else None)
let fG n i e=let g=n+i in let mutable n,l,p=n,1L,1L
(fun()->n<-n*G.[e]; p<-p*G.[e]; l<-l+p; let i=g*l-n in if i<=(int64 a) then Some(n,i,e) else None)
let rec fL n g=match n() with Some(f,i,e)->N.[(int i)]<-false; fL n ((fN f i (e+1))::g)
|_->match g with n::t->match n() with Some (n,i,e)->fL (fG n i e) g |_->fL n t
|_->N.[0]<-false; N
fL (fG 1L 0L 0) [fN 1L 0L 1]
 

The Task[edit]

Less than 2000
 
uT 2000|>Array.mapi(fun n g->(n,g))|>Array.filter(fun(_,n)->n)|>Array.chunkBySize 30|>Array.iter(fun n->n|>Array.iter(fst>>printf "%5d");printfn "")
 
Output:
    2    5   52   88   96  120  124  146  162  188  206  210  216  238  246  248  262  268  276  288  290  292  304  306  322  324  326  336  342  372
  406  408  426  430  448  472  474  498  516  518  520  530  540  552  556  562  576  584  612  624  626  628  658  668  670  708  714  718  726  732
  738  748  750  756  766  768  782  784  792  802  804  818  836  848  852  872  892  894  896  898  902  926  934  936  964  966  976  982  996 1002
 1028 1044 1046 1060 1068 1074 1078 1080 1102 1116 1128 1134 1146 1148 1150 1160 1162 1168 1180 1186 1192 1200 1212 1222 1236 1246 1248 1254 1256 1258
 1266 1272 1288 1296 1312 1314 1316 1318 1326 1332 1342 1346 1348 1360 1380 1388 1398 1404 1406 1418 1420 1422 1438 1476 1506 1508 1510 1522 1528 1538
 1542 1566 1578 1588 1596 1632 1642 1650 1680 1682 1692 1716 1718 1728 1732 1746 1758 1766 1774 1776 1806 1816 1820 1822 1830 1838 1840 1842 1844 1852
 1860 1866 1884 1888 1894 1896 1920 1922 1944 1956 1958 1960 1962 1972 1986 1992
Count less than or equal 100000
 
printfn "%d" (uT 100000|>Array.filter id|>Array.length)
 
Output:
13863
Real: 00:00:02.784, CPU: 00:00:02.750
Count less than or equal 1000000
 
printfn "%d" (uT 1000000|>Array.filter id|>Array.length)
 
Output:
150232
Real: 00:03:08.929, CPU: 00:03:08.859
Count less than or equal 2000000
 
printfn "%d" (uT 2000000|>Array.filter id|>Array.length)
 
Output:
305290
Real: 00:11:22.325, CPU: 00:11:21.828
Count less than or equal 3000000
 
<lang fsharp>
printfn "%d" (uT 3000000|>Array.filter id|>Array.length)
 
Output:
462110
Real: 00:24:00.126, CPU: 00:23:59.203

Go[edit]

This was originally based on the Wren example but has been modified somewhat to find untouchable numbers up to 1 million rather than 100,000 in a reasonable time. On my machine, the former (with a sieve factor of 63) took 31 minutes 9 seconds and the latter (with a sieve factor of 14) took 6.2 seconds.

package main
 
import "fmt"
 
func sumDivisors(n int) int {
sum := 1
k := 2
if n%2 == 0 {
k = 1
}
for i := 1 + k; i*i <= n; i += k {
if n%i == 0 {
sum += i
j := n / i
if j != i {
sum += j
}
}
}
return sum
}
 
func sieve(n int) []bool {
n++
s := make([]bool, n+1) // all false by default
for i := 6; i <= n; i++ {
sd := sumDivisors(i)
if sd <= n {
s[sd] = true
}
}
return s
}
 
func primeSieve(limit int) []bool {
limit++
// True denotes composite, false denotes prime.
c := make([]bool, limit) // all false by default
c[0] = true
c[1] = true
// no need to bother with even numbers over 2 for this task
p := 3 // Start from 3.
for {
p2 := p * p
if p2 >= limit {
break
}
for i := p2; i < limit; i += 2 * p {
c[i] = true
}
for {
p += 2
if !c[p] {
break
}
}
}
return c
}
 
func commatize(n int) string {
s := fmt.Sprintf("%d", n)
if n < 0 {
s = s[1:]
}
le := len(s)
for i := le - 3; i >= 1; i -= 3 {
s = s[0:i] + "," + s[i:]
}
if n >= 0 {
return s
}
return "-" + s
}
 
func main() {
limit := 1000000
c := primeSieve(limit)
s := sieve(63 * limit)
untouchable := []int{2, 5}
for n := 6; n <= limit; n += 2 {
if !s[n] && c[n-1] && c[n-3] {
untouchable = append(untouchable, n)
}
}
fmt.Println("List of untouchable numbers <= 2,000:")
count := 0
for i := 0; untouchable[i] <= 2000; i++ {
fmt.Printf("%6s", commatize(untouchable[i]))
if (i+1)%10 == 0 {
fmt.Println()
}
count++
}
fmt.Printf("\n\n%7s untouchable numbers were found <= 2,000\n", commatize(count))
p := 10
count = 0
for _, n := range untouchable {
count++
if n > p {
cc := commatize(count - 1)
cp := commatize(p)
fmt.Printf("%7s untouchable numbers were found <= %9s\n", cc, cp)
p = p * 10
if p == limit {
break
}
}
}
cu := commatize(len(untouchable))
cl := commatize(limit)
fmt.Printf("%7s untouchable numbers were found <= %s\n", cu, cl)
}
Output:
List of untouchable numbers <= 2,000:
     2     5    52    88    96   120   124   146   162   188
   206   210   216   238   246   248   262   268   276   288
   290   292   304   306   322   324   326   336   342   372
   406   408   426   430   448   472   474   498   516   518
   520   530   540   552   556   562   576   584   612   624
   626   628   658   668   670   708   714   718   726   732
   738   748   750   756   766   768   782   784   792   802
   804   818   836   848   852   872   892   894   896   898
   902   926   934   936   964   966   976   982   996 1,002
 1,028 1,044 1,046 1,060 1,068 1,074 1,078 1,080 1,102 1,116
 1,128 1,134 1,146 1,148 1,150 1,160 1,162 1,168 1,180 1,186
 1,192 1,200 1,212 1,222 1,236 1,246 1,248 1,254 1,256 1,258
 1,266 1,272 1,288 1,296 1,312 1,314 1,316 1,318 1,326 1,332
 1,342 1,346 1,348 1,360 1,380 1,388 1,398 1,404 1,406 1,418
 1,420 1,422 1,438 1,476 1,506 1,508 1,510 1,522 1,528 1,538
 1,542 1,566 1,578 1,588 1,596 1,632 1,642 1,650 1,680 1,682
 1,692 1,716 1,718 1,728 1,732 1,746 1,758 1,766 1,774 1,776
 1,806 1,816 1,820 1,822 1,830 1,838 1,840 1,842 1,844 1,852
 1,860 1,866 1,884 1,888 1,894 1,896 1,920 1,922 1,944 1,956
 1,958 1,960 1,962 1,972 1,986 1,992

    196 untouchable numbers were found  <=     2,000
      2 untouchable numbers were found  <=        10
      5 untouchable numbers were found  <=       100
     89 untouchable numbers were found  <=     1,000
  1,212 untouchable numbers were found  <=    10,000
 13,863 untouchable numbers were found  <=   100,000
150,232 untouchable numbers were found  <= 1,000,000

J[edit]

 
factor=: 3 : 0 NB. explicit
'primes powers'=. __&q: y
input_to_cartesian_product=. primes ^&.> i.&.> >: powers
cartesian_product=. , { input_to_cartesian_product
, */&> cartesian_product
)
 
factor=: [: , [: */&> [: { [: (^&.> i.&.>@>:)/ __&q: NB. tacit
 
 
proper_divisors=: [: }: factor
sum_of_proper_divisors=: +/@proper_divisors
 
candidates=: 5 , [: +: [: #\@i. >[email protected]: NB. within considered range, all but one candidate are even.
spds=:([:sum_of_proper_divisors"0(#\@i.-.i.&.:(p:inv))@*:)f. NB. remove primes which contribute 1
 

We may revisit to correct the larger untouchable tallies. The straightforward algorithm presented is a little bit efficient, and, I claim, the verb (candidates-.spds) produces the full list of untouchable numbers up to its argument. It considers the sum of proper divisors through the argument squared, less primes. Since J is an algorithm description language, it may be "fairer" to state in J that "more resources required" than in some other language. [1]

Time (seconds) and space (bytes) on a high end six year old (new in 2015) laptop computer.

   timespacex'([email protected]#)/:~~.spds 10000'
600.285 4.29497e9
   UNTOUCHABLES=:(candidates-.spds) 2000        NB. compute untouchable numbers

   /:~factor#UNTOUCHABLES                       NB. look for nice display size
1 2 4 7 14 28 49 98 196

   _14[\UNTOUCHABLES
   5    2   52   88   96  120  124  146  162  188  206  210  216  238
 246  248  262  268  276  288  290  292  304  306  322  324  326  336
 342  372  406  408  426  430  448  472  474  498  516  518  520  530
 540  552  556  562  576  584  612  624  626  628  658  668  670  708
 714  718  726  732  738  748  750  756  766  768  782  784  792  802
 804  818  836  848  852  872  892  894  896  898  902  926  934  936
 964  966  976  982  996 1002 1028 1044 1046 1060 1068 1074 1078 1080
1102 1116 1128 1134 1146 1148 1150 1160 1162 1168 1180 1186 1192 1200
1212 1222 1236 1246 1248 1254 1256 1258 1266 1272 1288 1296 1312 1314
1316 1318 1326 1332 1342 1346 1348 1360 1380 1388 1398 1404 1406 1418
1420 1422 1438 1476 1506 1508 1510 1522 1528 1538 1542 1566 1578 1588
1596 1632 1642 1650 1680 1682 1692 1716 1718 1728 1732 1746 1758 1766
1774 1776 1806 1816 1820 1822 1830 1838 1840 1842 1844 1852 1860 1866
1884 1888 1894 1896 1920 1922 1944 1956 1958 1960 1962 1972 1986 1992

   SBD=:spds 10000  NB. sums of proper divisors of many integers
   T=:/:~~.SBD  NB. sort the nub
   #T           NB. leaving this many touchable numbers, which are dense through at least 10000
47269787

   U=:([email protected]#)T NB. possible untouchable numbers

   head=: 'n';'tally of possible untouchable numbers to n'
   head,:,.L:0;/|:(,.U&I.)10^#\i.5  NB. testing to ten thousand squared
┌──────┬──────────────────────────────────────────┐
│n     │tally of possible untouchable numbers to n│
├──────┼──────────────────────────────────────────┤
│    10│    2                                     │
│   100│    5                                     │
│  1000│   89                                     │
│ 10000│ 1212                                     │
│100000│17538                                     │
└──────┴──────────────────────────────────────────┘

Julia[edit]

I can prove that the number to required to sieve to assure only untouchables for the interval 1:N is less than (N/2 - 1)^2 for larger N, but the 512,000,000 sieved below is just from doubling 1,000,000 and running the sieve until we get 150232 for the number of untouchables under 1,000,000.

using Primes
 
function properfactorsum(n)
f = [one(n)]
for (p,e) in factor(n)
f = reduce(vcat, [f*p^j for j in 1:e], init=f)
end
pop!(f)
return sum(f)
end
 
const maxtarget, sievelimit = 1_000_000, 512_000_000
const untouchables = ones(Bool, maxtarget)
 
for i in 2:sievelimit
n = properfactorsum(i)
if n <= maxtarget
untouchables[n] = false
end
end
for i in 6:maxtarget
if untouchables[i] && (isprime(i - 1) || isprime(i - 3))
untouchables[i] = false
end
end
 
println("The untouchable numbers ≤ 2000 are: ")
for (i, n) in enumerate(filter(x -> untouchables[x], 1:2000))
print(rpad(n, 5), i % 10 == 0 || i == 196 ? "\n" : "")
end
for N in [2000, 10, 100, 1000, 10_000, 100_000, 1_000_000]
println("The count of untouchable numbers ≤ $N is: ", count(x -> untouchables[x], 1:N))
end
 
Output:
The untouchable numbers ≤ 2000 are:
2    5    52   88   96   120  124  146  162  188
206  210  216  238  246  248  262  268  276  288
290  292  304  306  322  324  326  336  342  372
406  408  426  430  448  472  474  498  516  518
520  530  540  552  556  562  576  584  612  624
626  628  658  668  670  708  714  718  726  732
738  748  750  756  766  768  782  784  792  802
804  818  836  848  852  872  892  894  896  898
902  926  934  936  964  966  976  982  996  1002
1028 1044 1046 1060 1068 1074 1078 1080 1102 1116
1128 1134 1146 1148 1150 1160 1162 1168 1180 1186
1192 1200 1212 1222 1236 1246 1248 1254 1256 1258
1266 1272 1288 1296 1312 1314 1316 1318 1326 1332
1342 1346 1348 1360 1380 1388 1398 1404 1406 1418
1420 1422 1438 1476 1506 1508 1510 1522 1528 1538
1542 1566 1578 1588 1596 1632 1642 1650 1680 1682
1692 1716 1718 1728 1732 1746 1758 1766 1774 1776
1806 1816 1820 1822 1830 1838 1840 1842 1844 1852
1860 1866 1884 1888 1894 1896 1920 1922 1944 1956
1958 1960 1962 1972 1986 1992
The count of untouchable numbers ≤ 2000 is: 196
The count of untouchable numbers ≤ 10 is: 2
The count of untouchable numbers ≤ 100 is: 5
The count of untouchable numbers ≤ 1000 is: 89
The count of untouchable numbers ≤ 10000 is: 1212
The count of untouchable numbers ≤ 100000 is: 13863
The count of untouchable numbers ≤ 1000000 is: 150232

Perl[edit]

Library: ntheory
use strict;
use warnings;
use enum qw(False True);
use ntheory qw/divisor_sum is_prime/;
 
sub sieve {
my($n) = @_;
my %s;
for my $k (0 .. $n+1) {
my $sum = divisor_sum($k) - $k;
$s{$sum} = True if $sum <= $n+1;
}
%s
}
 
my(%s,%c);
my($max, $limit, $cnt) = (2000, 1e5, 0);
 
%s = sieve 14 * $limit;
!is_prime($_) and $c{$_} = True for 1..$limit;
my @untouchable = (2, 5);
for ( my $n = 6; $n <= $limit; $n += 2 ) {
push @untouchable, $n if !$s{$n} and $c{$n-1} and $c{$n-3};
}
map { $cnt++ if $_ <= $max } @untouchable;
print "Number of untouchable numbers ≤ $max : $cnt \n\n" .
(sprintf "@{['%6d' x $cnt]}", @untouchable[0..$cnt-1]) =~ s/(.{84})/$1\n/gr . "\n";
 
my($p, $count) = (10, 0);
my $fmt = "%6d untouchable numbers were found ≤ %7d\n";
for my $n (@untouchable) {
$count++;
if ($n > $p) {
printf $fmt, $count-1, $p;
printf($fmt, scalar @untouchable, $limit) and last if $limit == ($p *= 10)
}
}
Output:
Number of untouchable numbers ≤ 2000 : 196

     2     5    52    88    96   120   124   146   162   188   206   210   216   238
   246   248   262   268   276   288   290   292   304   306   322   324   326   336
   342   372   406   408   426   430   448   472   474   498   516   518   520   530
   540   552   556   562   576   584   612   624   626   628   658   668   670   708
   714   718   726   732   738   748   750   756   766   768   782   784   792   802
   804   818   836   848   852   872   892   894   896   898   902   926   934   936
   964   966   976   982   996  1002  1028  1044  1046  1060  1068  1074  1078  1080
  1102  1116  1128  1134  1146  1148  1150  1160  1162  1168  1180  1186  1192  1200
  1212  1222  1236  1246  1248  1254  1256  1258  1266  1272  1288  1296  1312  1314
  1316  1318  1326  1332  1342  1346  1348  1360  1380  1388  1398  1404  1406  1418
  1420  1422  1438  1476  1506  1508  1510  1522  1528  1538  1542  1566  1578  1588
  1596  1632  1642  1650  1680  1682  1692  1716  1718  1728  1732  1746  1758  1766
  1774  1776  1806  1816  1820  1822  1830  1838  1840  1842  1844  1852  1860  1866
  1884  1888  1894  1896  1920  1922  1944  1956  1958  1960  1962  1972  1986  1992

     2 untouchable numbers were found  ≤      10
     5 untouchable numbers were found  ≤     100
    89 untouchable numbers were found  ≤    1000
  1212 untouchable numbers were found  ≤   10000
 13863 untouchable numbers were found  ≤  100000

Phix[edit]

constant limz = {1,1,8,9,18,64} -- found by experiment
 
procedure untouchable(integer n, cols=0, tens=0)
atom t0 = time(), t1 = t0+1
bool tell = n>0
n = abs(n)
sequence sums = repeat(0,n+3)
for i=1 to n do
integer p = get_prime(i)
if p>n then exit end if
sums[p+1] = 1
sums[p+3] = 1
end for
sums[5] = 0
 
integer m = floor(log10(n))
integer lim = limz[m]*n
for j=2 to lim do
integer y = sum(factors(j,-1))
if y<=n then
sums[y] = 1
end if
if time()>t1 then
progress("j:%,d/%,d (%3.2f%%)\r",{j,lim,(j/lim)*100})
t1 = time()+1
end if
end for
progress("")
if tell then
printf(1,"The list of all untouchable numbers <= %d:\n",{n})
end if
string line = " 2 5"
integer cnt = 2
for t=6 to n by 2 do
if sums[t]=0 then
cnt += 1
if tell then
line &= sprintf("%,8d",t)
if remainder(cnt,cols)=0 then
printf(1,"%s\n",line)
line = ""
end if
end if
end if
end for
if tell then
if line!="" then
printf(1,"%s\n",line)
end if
printf(1,"\n")
end if
t0 = time()-t0
string t = iff(t0>1?elapsed(t0,fmt:=" (%s)"),"")
printf(1,"%,20d untouchable numbers were found <= %,d%s\n",{cnt,n,t})
for p=1 to tens do
untouchable(-power(10,p))
end for
end procedure
 
untouchable(2000, 10, 6)
Output:
The list of all untouchable numbers <= 2000:
       2       5      52      88      96     120     124     146     162     188
     206     210     216     238     246     248     262     268     276     288
     290     292     304     306     322     324     326     336     342     372
     406     408     426     430     448     472     474     498     516     518
     520     530     540     552     556     562     576     584     612     624
     626     628     658     668     670     708     714     718     726     732
     738     748     750     756     766     768     782     784     792     802
     804     818     836     848     852     872     892     894     896     898
     902     926     934     936     964     966     976     982     996   1,002
   1,028   1,044   1,046   1,060   1,068   1,074   1,078   1,080   1,102   1,116
   1,128   1,134   1,146   1,148   1,150   1,160   1,162   1,168   1,180   1,186
   1,192   1,200   1,212   1,222   1,236   1,246   1,248   1,254   1,256   1,258
   1,266   1,272   1,288   1,296   1,312   1,314   1,316   1,318   1,326   1,332
   1,342   1,346   1,348   1,360   1,380   1,388   1,398   1,404   1,406   1,418
   1,420   1,422   1,438   1,476   1,506   1,508   1,510   1,522   1,528   1,538
   1,542   1,566   1,578   1,588   1,596   1,632   1,642   1,650   1,680   1,682
   1,692   1,716   1,718   1,728   1,732   1,746   1,758   1,766   1,774   1,776
   1,806   1,816   1,820   1,822   1,830   1,838   1,840   1,842   1,844   1,852
   1,860   1,866   1,884   1,888   1,894   1,896   1,920   1,922   1,944   1,956
   1,958   1,960   1,962   1,972   1,986   1,992

                 196 untouchable numbers were found <= 2,000
                   2 untouchable numbers were found <= 10
                   6 untouchable numbers were found <= 100
                  89 untouchable numbers were found <= 1,000
               1,212 untouchable numbers were found <= 10,000
              13,863 untouchable numbers were found <= 100,000 (12.4s)
             150,232 untouchable numbers were found <= 1,000,000 (33 minutes and 55s)

for comparison, on the same box, the Julia entry took 58 minutes and 40s.

Raku[edit]

Borrow the proper divisors routine from here.

Translation of: Wren
# 20210220 Raku programming solution
 
sub propdiv (\x) {
my @l = 1 if x > 1;
(2 .. x.sqrt.floor).map: -> \d {
unless x % d { @l.push: d; my \y = x div d; @l.push: y if y != d }
}
@l
}
 
sub sieve (\n) {
my %s;
for (0..(n+1)) -> \k {
given ( [+] propdiv k ) { %s{$_} = True if $_(n+1) }
}
%s;
}
 
my \limit = 1e5;
my %c = ( grep { !.is-prime }, 1..limit ).Set; # store composites
my %s = sieve(14 * limit);
my @untouchable = 2, 5;
 
loop ( my \n = $ = 6 ; n ≤ limit ; n += 2 ) {
@untouchable.append(n) if (!%s{n} && %c{n-1} && %c{n-3})
}
 
my ($c, $last) = 0, False;
for @untouchable.rotor(10) {
say [~] @_».&{$c++ ; $_ > 2000 ?? ( $last = True and last ) !! .fmt: "%6d "}
$c-- and last if $last
}
 
say "\nList of untouchable numbers ≤ 2,000 : $c \n";
 
my ($p, $count) = 10,0;
BREAK: for @untouchable -> \n {
$count++;
if (n > $p) {
printf "%6d untouchable numbers were found ≤ %7d\n", $count-1, $p;
last BREAK if limit == ($p *= 10)
}
}
printf "%6d untouchable numbers were found ≤ %7d\n", +@untouchable, limit
Output:
     2      5     52     88     96    120    124    146    162    188
   206    210    216    238    246    248    262    268    276    288
   290    292    304    306    322    324    326    336    342    372
   406    408    426    430    448    472    474    498    516    518
   520    530    540    552    556    562    576    584    612    624
   626    628    658    668    670    708    714    718    726    732
   738    748    750    756    766    768    782    784    792    802
   804    818    836    848    852    872    892    894    896    898
   902    926    934    936    964    966    976    982    996   1002
  1028   1044   1046   1060   1068   1074   1078   1080   1102   1116
  1128   1134   1146   1148   1150   1160   1162   1168   1180   1186
  1192   1200   1212   1222   1236   1246   1248   1254   1256   1258
  1266   1272   1288   1296   1312   1314   1316   1318   1326   1332
  1342   1346   1348   1360   1380   1388   1398   1404   1406   1418
  1420   1422   1438   1476   1506   1508   1510   1522   1528   1538
  1542   1566   1578   1588   1596   1632   1642   1650   1680   1682
  1692   1716   1718   1728   1732   1746   1758   1766   1774   1776
  1806   1816   1820   1822   1830   1838   1840   1842   1844   1852
  1860   1866   1884   1888   1894   1896   1920   1922   1944   1956
  1958   1960   1962   1972   1986   1992

List of untouchable numbers ≤ 2,000 : 196

     2 untouchable numbers were found  ≤      10
     5 untouchable numbers were found  ≤     100
    89 untouchable numbers were found  ≤    1000
  1212 untouchable numbers were found  ≤   10000
 13863 untouchable numbers were found  ≤  100000

REXX[edit]

Some optimization was done to the code to produce prime numbers,   since a simple test could be made to exclude
the calculation of touchability for primes as the aliquot sum of a prime is always unity.
This saved around   15%   of the running time.

A fair amount of code was put into the generation of primes,   but the computation of the aliquot sum is the area
that consumes the most CPU time.

The REXX code below will accurately calculate   untouchable numbers   up to and including   100,000.   Beyond that,
a higher overreach   (the over option)   would be need to specified.

This version of REXX would need a 64-bit version to calculate the number of untouchable numbers for one million.

/*REXX pgm finds N untouchable numbers (numbers that can't be equal to any aliquot sum).*/
parse arg n cols tens over . /*obtain optional arguments from the CL*/
if n='' | n=="," then n=2000 /*Not specified? Then use the default.*/
if cols='' | cols=="," | cols==0 then cols= 10 /* " " " " " " */
if tens='' | tens=="," then tens= 0 /* " " " " " " */
if over='' | over=="," then over= 20 /* " " " " " " */
tell= n>0; n= abs(n) /*N>0? Then display the untouchable #s*/
call genP n * over /*call routine to generate some primes.*/
u.= 0 /*define all possible aliquot sums ≡ 0.*/
do p=1 for #; _= @.p + 1; u._= 1 /*any prime+1 is not an untouchable.*/
_= @.p + 3; u._= 1 /* " prime+3 " " " " */
end /*p*/ /* [↑] this will also rule out 5. */
u.5= 0 /*special case as prime 2 + 3 sum to 5.*/
do j=2 for lim; if !.j then iterate /*Is J a prime? Yes, then skip it. */
y= sigmaP() /*compute: aliquot sum (sigma P) of J.*/
if y<=n then u.y= 1 /*mark Y as a touchable if in range. */
end /*j*/
call show /*maybe show untouchable #s and a count*/
if tens>0 then call powers /*Any "tens" specified? Calculate 'em.*/
exit cnt /*stick a fork in it, we're all done. */
/*──────────────────────────────────────────────────────────────────────────────────────*/
commas: parse arg ?; do jc=length(?)-3 to 1 by -3; ?=insert(',', ?, jc); end; return ?
genSq: do _=1 until _*_>lim; q._= _*_; end; q._= _*_; _= _+1; q._= _*_; return
grid: $= $ right( commas(t), w); if cnt//cols==0 then do; say $; $=; end; return
powers: do pr=1 for tens; call 'UNTOUCHA' -(10**pr); end /*recurse*/; return
/*──────────────────────────────────────────────────────────────────────────────────────*/
genP: #= 9; @.1=2; @.2=3; @.3=5; @.4=7; @.5=11; @.6=13; @.7=17; @.8=19; @.9=23 /*a list*/
 !.=0;  !.2=1; !.3=1; !.5=1; !.7=1; !.11=1; !.13=1; !.17=1; !.19=1  !.23=1 /*primes*/
parse arg lim; call genSq /*define the (high) limit for searching*/
qq.10= 100 /*define square of the 10th prime index*/
do [email protected].#+6 by 2 to lim /*find odd primes from here on forward.*/
parse var j '' -1 _; if _==5 then iterate; if j// 3==0 then iterate
if j// 7==0 then iterate; if j//11==0 then iterate; if j//13==0 then iterate
if j//17==0 then iterate; if j//19==0 then iterate; if j//23==0 then iterate
/*start dividing by the tenth prime: 29*/
do k=10 while qq.k <= j /* [↓] divide J by known odd primes.*/
if j//@.k==0 then iterate j /*J ÷ by a prime? Then ¬prime. ___ */
end /*k*/ /* [↑] only process numbers ≤ √ J */
#= #+1; @.#= j /*bump prime count; assign a new prime.*/
 !.j= 1; qq.#= j*j /*mark prime; compute square of prime.*/
end /*j*/; return /*#: is the number of primes generated*/
/*──────────────────────────────────────────────────────────────────────────────────────*/
show: w=7; $= right(2, w+1) right(5, w) /*start the list of an even prime and 5*/
cnt= 2 /*count of the only two primes in list.*/
do t=6 by 2 to n; if u.t then iterate /*Is T touchable? Then skip it. */
cnt= cnt + 1; if tell then call grid /*bump count; maybe show a grid line. */
end /*t*/
if tell & $\=='' then say $ /*display a residual grid line, if any.*/
if tell then say /*show a spacing blank line for output.*/
if n>0 then say right( commas(cnt), 20) , /*indent the output a bit.*/
' untouchable numbers were found ≤ ' commas(n); return
/*──────────────────────────────────────────────────────────────────────────────────────*/
sigmaP: s= 1 /*set initial sigma sum (S) to 1. ___*/
if j//2 then do m=3 by 2 while q.m<j /*divide by odd integers up to the √ J */
if j//m==0 then s=s+m+j%m /*add the two divisors to the sum. */
end /*m*/ /* [↑] process an odd integer. ___*/
else do m=2 while q.m<j /*divide by all integers up to the √ J */
if j//m==0 then s=s+m+j%m /*add the two divisors to the sum. */
end /*m*/ /* [↑] process an even integer. ___*/
if q.m==j then return s + m /*Was J a square? If so, add √ J */
return s /* No, just return. */
output   when using the default inputs:
       2       5      52      88      96     120     124     146     162     188
     206     210     216     238     246     248     262     268     276     288
     290     292     304     306     322     324     326     336     342     372
     406     408     426     430     448     472     474     498     516     518
     520     530     540     552     556     562     576     584     612     624
     626     628     658     668     670     708     714     718     726     732
     738     748     750     756     766     768     782     784     792     802
     804     818     836     848     852     872     892     894     896     898
     902     926     934     936     964     966     976     982     996   1,002
   1,028   1,044   1,046   1,060   1,068   1,074   1,078   1,080   1,102   1,116
   1,128   1,134   1,146   1,148   1,150   1,160   1,162   1,168   1,180   1,186
   1,192   1,200   1,212   1,222   1,236   1,246   1,248   1,254   1,256   1,258
   1,266   1,272   1,288   1,296   1,312   1,314   1,316   1,318   1,326   1,332
   1,342   1,346   1,348   1,360   1,380   1,388   1,398   1,404   1,406   1,418
   1,420   1,422   1,438   1,476   1,506   1,508   1,510   1,522   1,528   1,538
   1,542   1,566   1,578   1,588   1,596   1,632   1,642   1,650   1,680   1,682
   1,692   1,716   1,718   1,728   1,732   1,746   1,758   1,766   1,774   1,776
   1,806   1,816   1,820   1,822   1,830   1,838   1,840   1,842   1,844   1,852
   1,860   1,866   1,884   1,888   1,894   1,896   1,920   1,922   1,944   1,956
   1,958   1,960   1,962   1,972   1,986   1,992

                 196  untouchable numbers were found  ≤  2,000 
output   when using the inputs:     0   ,   5
                   2  untouchable numbers were found  ≤  10
                   5  untouchable numbers were found  ≤  100
                  89  untouchable numbers were found  ≤  1,000
               1,212  untouchable numbers were found  ≤  10,000
              13,863  untouchable numbers were found  ≤  100,000

Wren[edit]

Library: Wren-seq
Library: Wren-math
Library: Wren-fmt

Not an easy task as, even allowing for the prime tests, it's difficult to know how far you need to sieve to get the right answers. The parameters here were found by trial and error.

import "/math" for Int, Nums
import "/seq" for Lst
import "/fmt" for Fmt
 
var sieve = Fn.new { |n|
n = n + 1
var s = List.filled(n+1, false)
for (i in 0..n) {
var sum = Nums.sum(Int.properDivisors(i))
if (sum <= n) s[sum] = true
}
return s
}
 
var limit = 1e5
var c = Int.primeSieve(limit, false)
var s = sieve.call(14 * limit)
var untouchable = [2, 5]
var n = 6
while (n <= limit) {
if (!s[n] && c[n-1] && c[n-3]) untouchable.add(n)
n = n + 2
}
 
System.print("List of untouchable numbers <= 2,000:")
for (chunk in Lst.chunks(untouchable.where { |n| n <= 2000 }.toList, 10)) {
Fmt.print("$,6d", chunk)
}
System.print()
Fmt.print("$,6d untouchable numbers were found <= 2,000", untouchable.count { |n| n <= 2000 })
var p = 10
var count = 0
for (n in untouchable) {
count = count + 1
if (n > p) {
Fmt.print("$,6d untouchable numbers were found <= $,7d", count-1, p)
p = p * 10
if (p == limit) break
}
}
Fmt.print("$,6d untouchable numbers were found <= $,d", untouchable.count, limit)
Output:
List of untouchable numbers <= 2,000:
     2      5     52     88     96    120    124    146    162    188
   206    210    216    238    246    248    262    268    276    288
   290    292    304    306    322    324    326    336    342    372
   406    408    426    430    448    472    474    498    516    518
   520    530    540    552    556    562    576    584    612    624
   626    628    658    668    670    708    714    718    726    732
   738    748    750    756    766    768    782    784    792    802
   804    818    836    848    852    872    892    894    896    898
   902    926    934    936    964    966    976    982    996  1,002
 1,028  1,044  1,046  1,060  1,068  1,074  1,078  1,080  1,102  1,116
 1,128  1,134  1,146  1,148  1,150  1,160  1,162  1,168  1,180  1,186
 1,192  1,200  1,212  1,222  1,236  1,246  1,248  1,254  1,256  1,258
 1,266  1,272  1,288  1,296  1,312  1,314  1,316  1,318  1,326  1,332
 1,342  1,346  1,348  1,360  1,380  1,388  1,398  1,404  1,406  1,418
 1,420  1,422  1,438  1,476  1,506  1,508  1,510  1,522  1,528  1,538
 1,542  1,566  1,578  1,588  1,596  1,632  1,642  1,650  1,680  1,682
 1,692  1,716  1,718  1,728  1,732  1,746  1,758  1,766  1,774  1,776
 1,806  1,816  1,820  1,822  1,830  1,838  1,840  1,842  1,844  1,852
 1,860  1,866  1,884  1,888  1,894  1,896  1,920  1,922  1,944  1,956
 1,958  1,960  1,962  1,972  1,986  1,992

   196 untouchable numbers were found  <=   2,000
     2 untouchable numbers were found  <=      10
     5 untouchable numbers were found  <=     100
    89 untouchable numbers were found  <=   1,000
 1,212 untouchable numbers were found  <=  10,000
13,863 untouchable numbers were found  <= 100,000