Untouchable numbers

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   p²   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

•   Wolfram MathWorld:   untouchable number.
•   OEIS:   A005114 untouchable numbers.
•   OEIS:   a list of all untouchable numbers below 100,000   (inclusive).
•   Wikipedia: untouchable number.
•   Wikipedia: Goldbach's conjecture.

ALGOL 68

Generates the proper divisor sums with a sieve-like process.
As noted by the Wren, Go, etc. samples, it is hard to determine what upper limit to use to eliminate non-untouchable numbers. It seems that using the proper divisor sums only, to find untouchables up to limit, limit^2 must be considered.
However, using the observations of the Wren, Go etc. solutions that by using the facts prime + 1 and prime + 3 are not untouchable, and that (probably) 5 is the only odd untouchable, the untouchables up to 1 000 000 can be found by considering numbers up to 64 000 000 (possibly less could be used). At least, this sample finds the same values as the other ones :).
Possibly, this works because the prime + 1 value eliminates the need to calculate the proper divisor sum of p^2 which appears to help a lot when p is close to the upper limit. Why the 64 * limit (or 63 * limit) is valid is unclear (to me, anyway).

Note that under Windows (and possibly under Linux), Algol 68G requires that the heap size be increased in order to allow arrays big enough to handle 100 000 and 1 000 000 untouchable numbers. See ALGOL_68_Genie#Using_a_Large_Heap.

Library: ALGOL 68-primes

BEGIN # find some untouchable numbers - numbers not equal to the sum of the   #
      # proper divisors of any +ve integer                                    #
    INT max untouchable = 1 000 000;
    # a table of the untouchable numbers                                      #
    [ 1 : max untouchable ]BOOL untouchable; FOR i TO UPB untouchable DO untouchable[ i ] := TRUE OD;
    # show the counts of untouchable numbers found                            #
    PROC show untouchable statistics = VOID:
         BEGIN
            print( ( "Untouchable numbers:", newline ) );
            INT u count := 0;
            FOR i TO UPB untouchable DO
                IF untouchable[ i ] THEN u count +:= 1 FI;
                IF i =        10
                OR i =       100
                OR i =     1 000
                OR i =    10 000
                OR i =   100 000
                OR i = 1 000 000
                THEN
                    print( ( whole( u count, -7 ), " to ", whole( i, -8 ), newline ) )
                FI
            OD
         END; # show untouchable counts #
    # prints the untouchable numbers up to n                                  #
    PROC print untouchables = ( INT n )VOID:
         BEGIN
            print( ( "Untouchable numbers up to ", whole( n, 0 ), newline ) );
            INT u count := 0;
            FOR i TO n DO
                IF untouchable[ i ] THEN
                    print( ( whole( i, -4 ) ) );
                    IF u count +:= 1;
                       u count MOD 16 = 0
                    THEN print( ( newline ) )
                    ELSE print( ( " " ) )
                    FI
                FI
            OD;
            print( ( newline ) );
            print( ( whole( u count, -7 ), " to ", whole( n, -8 ), newline ) )
         END; # print untouchables #
    # find the untouchable numbers                                            #
    # to find untouchable numbers up to e.g.: 10 000, we need to sieve up to  #
    # 10 000 ^2 i.e. 100 000 000                                              #
    # however if we also use the facts that no untouchable = prime + 1        #
    # and no untouchable = odd prime + 3 and 5 is (very probably) the only    #
    # odd untouchable, other samples suggest we can use limit * 64 to find    #
    # untlouchables up to 1 000 000 - experimentation reveals this to be true #
    # assume the conjecture that there are no odd untouchables except 5       #
    BEGIN
        untouchable[ 1 ] := FALSE;
        untouchable[ 3 ] := FALSE;
        FOR i FROM 7 BY 2 TO UPB untouchable DO untouchable[ i ] := FALSE OD
    END;
    # sieve the primes to max untouchable and flag the non untouchables       #
    BEGIN
        PR read "primes.incl.a68" PR
        []BOOL prime = PRIMESIEVE max untouchable;
        FOR i FROM 3 BY 2 TO UPB prime DO
            IF prime[ i ] THEN
                IF i < max untouchable THEN
                    untouchable[ i + 1 ] := FALSE;
                    IF i < ( max untouchable - 2 ) THEN
                        untouchable[ i + 3 ] := FALSE
                    FI
                FI
            FI
        OD;
        untouchable[ 2 + 1 ] := FALSE # special case for the only even prime  #
    END;
    # construct the proper divisor sums and flag the non untouchables         #
    BEGIN
        [ 1 : max untouchable * 64 ]INT spd;
        FOR i TO UPB spd DO spd[ i ] := 1 OD;
        FOR i FROM 2 TO UPB spd DO
            FOR j FROM i + i BY i TO UPB spd DO spd[ j ] +:= i OD
        OD;
        FOR i TO UPB spd DO
            IF spd[ i ] <= UPB untouchable THEN untouchable[ spd[ i ] ] := FALSE FI
        OD
    END;
    # show the untouchable numbers up to 2000                                 #
    print untouchables( 2 000 );
    # show the counts of untouchable numbers                                  #
    show untouchable statistics
END

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
    196 to     2000
Untouchable numbers:
      2 to       10
      5 to      100
     89 to     1000
   1212 to    10000
  13863 to   100000
 150232 to  1000000

C

Run time around 14 seconds on my Core i7 machine.

#include <stdio.h>
#include <stdlib.h>
#include <stdbool.h>
#include <locale.h>

bool *primeSieve(int limit) {
    int i, p;
    limit++;
    // True denotes composite, false denotes prime.
    bool *c = calloc(limit, sizeof(bool)); // all false by default
    c[0] = true;
    c[1] = true;
    for (i = 4; i < limit; i += 2) c[i] = true;
    p = 3; // Start from 3.
    while (true) {
        int p2 = p * p;
        if (p2 >= limit) break;
        for (i = p2; i < limit; i += 2 * p) c[i] = true;
        while (true) {
            p += 2;
            if (!c[p]) break;
        }
    }
    return c;
}

int main() {
    const int limit = 1000000;
    int i, j, n, uc = 2, p = 10, m = 63, ul = 151000;
    bool *c = primeSieve(limit);
    n = m * limit + 1;
    int *sumDivs = (int *)calloc(n, sizeof(int));
    for (i = 1; i < n; ++i) {
        for (j = i; j < n; j += i) sumDivs[j] += i;
    }
    bool *s = (bool *)calloc(n, sizeof(bool)); // all false
    for (i = 1; i < n; ++i) {
        int sum = sumDivs[i] - i; // proper divs sum
        if (sum <= n) s[sum] = true;
    } 
    free(sumDivs);
    int *untouchable = (int *)malloc(ul * sizeof(int));
    untouchable[0] = 2;
    untouchable[1] = 5;
    for (n = 6; n <= limit; n += 2) {
        if (!s[n] && c[n-1] && c[n-3]) untouchable[uc++] = n;
    }
    setlocale(LC_NUMERIC, "");
    printf("List of untouchable numbers <= 2,000:\n");
    for (i = 0; i < uc; ++i) {
        j = untouchable[i];
        if (j > 2000) break;
        printf("%'6d ", j);
        if (!((i+1) % 10)) printf("\n");
    }
    printf("\n\n%'7d untouchable numbers were found  <=     2,000\n", i);
    for (i = 0; i < uc; ++i) {
        j = untouchable[i];
        if (j > p) {
            printf("%'7d untouchable numbers were found  <= %'9d\n", i, p);
            p *= 10;
            if (p == limit) break;
        }
    }
    printf("%'7d untouchable numbers were found  <= %'d\n", uc, limit);
    free(c);
    free(s);
    free(untouchable);
    return 0;
}

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

C++

This solution implements Talk:Untouchable_numbers#Nice_recursive_solution

The Function

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

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

This example does not show the output mentioned in the task description on this page (or a page linked to from here). Please ensure that it meets all task requirements and remove this message. Note that phrases in task descriptions such as "print and display" and "print and show" for example, indicate that (reasonable length) output be a part of a language's solution.

Library: System.SysUtils

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.

FreeBASIC

Const max_untouchable = 1e6

Dim Shared untouchable(1 To max_untouchable) As Uinteger
For i As Uinteger = 1 To Ubound(untouchable)
    untouchable(i) = True
Next i

Sub show_untouchable_statistics()
    Dim As Uinteger i, cnt = 0
    For i = 1 To Ubound(untouchable)
        If untouchable(i) Then cnt += 1
        If i = 10 Orelse i = 100 Orelse i = 1000 Orelse i = 10000 Orelse i = 1e5 Then
            Print Using "###,### untouchable numbers were found  <= #,###,###"; cnt;  i
        End If
    Next i
End Sub

Sub print_untouchables(n As Uinteger)
    Print Using "List of untouchable numbers <= #,###:"; n
    Dim As Uinteger i, cnt = 0
    For i = 1 To n
        If untouchable(i) Then
            Print Using "##,###"; i;
            cnt += 1
            Print Iif(cnt Mod 16 = 0, !"\n", " ");
        End If
    Next i
    Print: Print
    Print Using "###,### untouchable numbers were found  <= #,###,###"; cnt; n
End Sub

Dim As Uinteger i, j
untouchable(1) = False
untouchable(3) = False
For i = 7 To Ubound(untouchable) Step 2
    untouchable(i) = False
Next i

Dim Shared spd(1 To max_untouchable * 64) As Uinteger
Dim As Uinteger ub = Ubound(spd)
For i = 1 To ub
    spd(i) = 1
Next i
For i = 2 To ub
    For j = i + i To ub Step i
        spd(j) += i
    Next j
Next i
For i = 1 To ub
    If spd(i) <= Ubound(untouchable) Then untouchable(spd(i)) = False
Next i

' Show the untouchable numbers up to 2000
print_untouchables(2000)
' Show the counts of untouchable numbers
show_untouchable_statistics()

Sleep

F#

The Function

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

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

printfn "%d" (uT 3000000|>Array.filter id|>Array.length)

Output:

462110
Real: 00:24:00.126, CPU: 00:23:59.203

Go

Version1

This was originally based on the Wren Version 1 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

Version 2

Translation of: C

Library: Go-rcu

This version uses a much more efficient algorithm for calculating the sums of divisors in bulk.

Run time for finding untouchable numbers up to 1 million is now only 11 seconds which (surprisingly) is faster than C itself on the same machine.

Also, since the first version was written, some of the functions used have now been incorporated into the above library.

package main

import (
    "fmt"
    "rcu"
)

func main() {
    limit := 1000000
    m := 63
    c := rcu.PrimeSieve(limit, false)
    n := m*limit + 1
    sumDivs := make([]int, n)
    for i := 1; i < n; i++ {
        for j := i; j < n; j += i {
            sumDivs[j] += i
        }
    }
    s := make([]bool, n) // all false
    for i := 1; i < n; i++ {
        sum := sumDivs[i] - i // proper divs sum
        if sum <= n {
            s[sum] = true
        }
    }
    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", rcu.Commatize(untouchable[i]))
        if (i+1)%10 == 0 {
            fmt.Println()
        }
        count++
    }
    fmt.Printf("\n\n%7s untouchable numbers were found  <=     2,000\n", rcu.Commatize(count))
    p := 10
    count = 0
    for _, n := range untouchable {
        count++
        if n > p {
            cc := rcu.Commatize(count - 1)
            cp := rcu.Commatize(p)
            fmt.Printf("%7s untouchable numbers were found  <= %9s\n", cc, cp)
            p = p * 10
            if p == limit {
                break
            }
        }
    }
    cu := rcu.Commatize(len(untouchable))
    cl := rcu.Commatize(limit)
    fmt.Printf("%7s untouchable numbers were found  <= %s\n", cu, cl)
}

Output:

Same as Version 1.

J

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. >.@-:  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'(-.candidates@#)/:~~.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=:(-.~candidates@#)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

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

Mathematica/Wolfram Language

f = DivisorSigma[1, #] - # &;
limit = 10^5;
c = Not /@ PrimeQ[Range[limit]];
slimit = 15 limit;
s = ConstantArray[False, slimit + 1];
untouchable = {2, 5};
Do[
 val = f[i];
 If[val <= slimit,
  s[[val]] = True
  ]
 ,
 {i, 6, slimit}
]
Do[
 If[! s[[n]],
  If[c[[n - 1]],
   If[c[[n - 3]],
    AppendTo[untouchable, n]
    ]
   ]
  ]
 ,
 {n, 6, limit, 2}
]
Multicolumn[Select[untouchable, LessEqualThan[2000]]]
Count[untouchable, _?(LessEqualThan[2000])]
Count[untouchable, _?(LessEqualThan[10])]
Count[untouchable, _?(LessEqualThan[100])]
Count[untouchable, _?(LessEqualThan[1000])]
Count[untouchable, _?(LessEqualThan[10000])]
Count[untouchable, _?(LessEqualThan[100000])]

Output:

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

Nim

I borrowed some ideas from Go version 1, but keep the limit to 100_000 as in Wren version 1.

import math, strutils

const
  Lim1 = 100_000    # Limit for untouchable numbers.
  Lim2 = 14 * Lim1  # Limit for computation of sum of divisors.

proc sumdiv(n: uint): uint =
  ## Return the sum of the strict divisors of "n".
  result = 1
  let r = sqrt(n.float).uint
  let k = if (n and 1) == 0: 1u else: 2u
  for d in countup(k + 1, r, k):
    if n mod d == 0:
      result += d
      let q = n div d
      if q != d: result += q

var
  isSumDiv: array[1..Lim2, bool]
  isPrime: array[1..Lim1, bool]

# Fill both sieves in a single pass.
for n in 1u..Lim2:
  let s = sumdiv(n)
  if s <= Lim2:
    isSumDiv[s] = true
    if s == 1 and n <= Lim1:
      isPrime[n] = true
isPrime[1] = false

# Build list of untouchable numbers.
var list = @[2, 5]
for n in countup(6, Lim1, 2):
  if not (isSumDiv[n] or isPrime[n - 1] or isPrime[n - 3]):
    list.add n

echo "Untouchable numbers ≤ 2000:"
var count, lcount = 0
for n in list:
  if n <= 2000:
    stdout.write ($n).align(5)
    inc count
    inc lcount
    if lcount == 20:
      echo()
      lcount = 0
  else:
    if lcount > 0: echo()
    break

const CountMessage = "There are $1 untouchable numbers ≤ $2."
echo CountMessage.format(count, 2000), '\n'

count = 0
var lim = 10
for n in list:
  if n > lim:
    echo CountMessage.format(count, lim)
    lim *= 10
  inc count
if lim == Lim1:
  # Emit last message.
  echo CountMessage.format(count, lim)

Output:

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 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
There are 196 untouchable numbers ≤ 2000.

There are 2 untouchable numbers ≤ 10.
There are 5 untouchable numbers ≤ 100.
There are 89 untouchable numbers ≤ 1000.
There are 1212 untouchable numbers ≤ 10000.
There are 13863 untouchable numbers ≤ 100000.

Pascal

modified Factors_of_an_integer#using_Prime_decomposition to calculate only sum of divisors
Appended a list of count of untouchable numbers out of math.dartmouth.edu/~carlp/uupaper3.pdf
Nigel is still right, that this procedure is not the way to get proven results.

program UntouchableNumbers;
program UntouchableNumbers;
{$IFDEF FPC}
  {$MODE DELPHI}  {$OPTIMIZATION ON,ALL}  {$COPERATORS ON}
  {$CODEALIGN proc=16,loop=4}
{$ELSE}
  {$APPTYPE CONSOLE}
{$ENDIF}
uses
  sysutils,strutils
{$IFDEF WINDOWS},Windows{$ENDIF}
  ;
const
  MAXPRIME = 1742537;
  //sqr(MaxPrime) = 3e12
  LIMIT =  5*1000*1000;
  LIMIT_mul = trunc(exp(ln(LIMIT)/3))+1;

const
  SizePrDeFe = 16*8192;//*size of(tprimeFac) =16 byte  2 Mb ~ level 3 cache
type
  tdigits = array [0..31] of Uint32;
  tprimeFac = packed record
                 pfSumOfDivs,
                 pfRemain : Uint64;
               end;
  tpPrimeFac = ^tprimeFac;

  tPrimeDecompField = array[0..SizePrDeFe-1] of tprimeFac;

  tPrimes = array[0..1 shl 17-1] of Uint32;

var
  {$ALIGN 16}
  PrimeDecompField :tPrimeDecompField;
  {$ALIGN 16}
  SmallPrimes: tPrimes;
  pdfIDX,pdfOfs: NativeInt;
  TD : Int64;

procedure OutCounts(pUntouch:pByte);
var
  n,cnt,lim,deltaLim : NativeInt;
Begin
  n := 0;
  cnt := 0;
  deltaLim := 100;
  lim := deltaLim;
  repeat
    cnt += 1-pUntouch[n];
    if n = lim then
    Begin
      writeln(Numb2USA(IntToStr(lim)):13,' ',Numb2USA(IntToStr(cnt)):12);
      lim += deltaLim;
      if lim = 10*deltaLim then
      begin
        deltaLim *=10;
        lim := deltaLim;
        writeln;
      end;
    end;

    inc(n);
  until n > LIMIT;
end;

function OutN(n:UInt64):UInt64;
begin
  write(Numb2USA(IntToStr(n)):15,' dt ',(GettickCount64-TD)/1000:5:3,' s'#13);
  TD := GettickCount64;
  result := n+LIMIT;
end;

//######################################################################
//gets sum of divisors of consecutive integers fast
procedure InitSmallPrimes;
//get primes. Sieving only odd numbers
var
  pr : array[0..MAXPRIME] of byte;
  p,j,d,flipflop :NativeUInt;
Begin
  SmallPrimes[0] := 2;
  fillchar(pr[0],SizeOf(pr),#0);
  p := 0;
  repeat
    repeat
      p +=1
    until pr[p]= 0;
    j := (p+1)*p*2;
    if j>MAXPRIME then
      BREAK;
    d := 2*p+1;
    repeat
      pr[j] := 1;
      j += d;
    until j>MAXPRIME;
  until false;

  SmallPrimes[1] := 3;
  SmallPrimes[2] := 5;
  j := 3;
  flipflop := (2+1)-1;//7+2*2->11+2*1->13 ,17 ,19 , 23
  p := 3;
  repeat
    if pr[p] = 0 then
    begin
      SmallPrimes[j] := 2*p+1;
      inc(j);
    end;
    p+=flipflop;
    flipflop := 3-flipflop;
  until (p > MAXPRIME) OR (j>High(SmallPrimes));
end;

function CnvtoBASE(var dgt:tDigits;n:Uint64;base:NativeUint):NativeInt;
//n must be multiple of base aka n mod base must be 0
var
  q,r: Uint64;
  i : NativeInt;
Begin
  fillchar(dgt,SizeOf(dgt),#0);
  i := 0;
  n := n div base;
  result := 0;
  repeat
    r := n;
    q := n div base;
    r  -= q*base;
    n := q;
    dgt[i] := r;
    inc(i);
  until (q = 0);
  //searching lowest pot in base
  result := 0;
  while (result<i) AND (dgt[result] = 0) do
    inc(result);
  inc(result);
end;

function IncByOneInBase(var dgt:tDigits;base:NativeInt):NativeInt;
var
  q :NativeInt;
Begin
  result := 0;
  q := dgt[result]+1;
  if q = base then
    repeat
      dgt[result] := 0;
      inc(result);
      q := dgt[result]+1;
    until q <> base;
  dgt[result] := q;
  result +=1;
end;

procedure CalcSumOfDivs(var pdf:tPrimeDecompField;var dgt:tDigits;n,k,pr:Uint64);
var
  fac,s :Uint64;
  j : Int32;
Begin
  //j is power of prime
  j := CnvtoBASE(dgt,n+k,pr);
  repeat
    fac := 1;
    s := 1;
    repeat
      fac *= pr;
      dec(j);
      s += fac;
    until j<= 0;
    with pdf[k] do
    Begin
      pfSumOfDivs *= s;
      pfRemain := pfRemain DIV fac;
    end;
    j := IncByOneInBase(dgt,pr);
    k += pr;
  until k >= SizePrDeFe;
end;

function SieveOneSieve(var pdf:tPrimeDecompField):boolean;
var
  dgt:tDigits;
  i,j,k,pr,n,MaxP : Uint64;
begin
  n := pdfOfs;
  if n+SizePrDeFe >= sqr(SmallPrimes[High(SmallPrimes)]) then
    EXIT(FALSE);
  //init
  for i := 0 to SizePrDeFe-1 do
  begin
    with pdf[i] do
    Begin
      pfSumOfDivs := 1;
      pfRemain := n+i;
    end;
  end;
  //first factor 2. Make n+i even
  i := (pdfIdx+n) AND 1;
  IF (n = 0) AND (pdfIdx<2)  then
    i := 2;

  repeat
    with pdf[i] do
    begin
      j := BsfQWord(n+i);
      pfRemain := (n+i) shr j;
      pfSumOfDivs := (Uint64(1) shl (j+1))-1;
    end;
    i += 2;
  until i >=SizePrDeFe;
  //i now index in SmallPrimes
  i := 0;
  maxP := trunc(sqrt(n+SizePrDeFe))+1;
  repeat
    //search next prime that is in bounds of sieve
    if n = 0 then
    begin
      repeat
        inc(i);
        pr := SmallPrimes[i];
        k := pr-n MOD pr;
        if k < SizePrDeFe then
          break;
      until pr > MaxP;
    end
    else
    begin
      repeat
        inc(i);
        pr := SmallPrimes[i];
        k := pr-n MOD pr;
        if (k = pr) AND (n>0) then
          k:= 0;
        if k < SizePrDeFe then
          break;
      until pr > MaxP;
    end;

    //no need to use higher primes
    if pr > maxP then
      BREAK;

    CalcSumOfDivs(pdf,dgt,n,k,pr);
  until false;

  //correct sum of & count of divisors
  for i := 0 to High(pdf) do
  Begin
    with pdf[i] do
    begin
      j := pfRemain;
      if j <> 1 then
        pfSumOFDivs *= (j+1);
    end;
  end;
  result := true;
end;

function NextSieve:boolean;
begin
  dec(pdfIDX,SizePrDeFe);
  inc(pdfOfs,SizePrDeFe);
  result := SieveOneSieve(PrimeDecompField);
end;

function GetNextPrimeDecomp:tpPrimeFac;
begin
  if pdfIDX >= SizePrDeFe then
    if Not(NextSieve) then
    Begin
      writeln('of limits ');
      EXIT(NIL);
    end;
  result := @PrimeDecompField[pdfIDX];
  inc(pdfIDX);
end;

function Init_Sieve(n:NativeUint):boolean;
//Init Sieve pdfIdx,pdfOfs are Global
begin
  pdfIdx := n MOD SizePrDeFe;
  pdfOfs := n-pdfIdx;
  result := SieveOneSieve(PrimeDecompField);
end;
//gets sum of divisors of consecutive integers fast
//######################################################################

procedure CheckRest(n: Uint64;pUntouch:pByte);
var
  k,lim : Uint64;
begin
  lim := 2*LIMIT;
  repeat
    k := GetNextPrimeDecomp^.pfSumOfDivs-n;
    inc(n);
    if Not(ODD(k)) AND (k<=LIMIT) then
      pUntouch[k] := 1;
  // showing still alive not for TIO.RUN
//    if n >= lim then  lim := OutN(n);
  until n >LIMIT_mul*LIMIT;
end;

var
  Untouch : array of byte;
  pUntouch: pByte;
  puQW  : pQword;
  T0:Int64;
  n,k : NativeInt;
Begin
  if sqrt(LIMIT_mul*LIMIT) >=MAXPRIME then
  Begin
    writeln('Need to extend count of primes > ',
      trunc(sqrt(LIMIT_mul*LIMIT))+1);
    HALT(0);
  end;

  setlength(untouch,LIMIT+8+1);
  pUntouch := @untouch[0];
  //Mark all odd as touchable
  puQW := @pUntouch[0];
  For n := 0 to LIMIT DIV 8 do puQW[n] := $0100010001000100;

  InitSmallPrimes;
  T0 := GetTickCount64;
  writeln('LIMIT = ',Numb2USA(IntToStr(LIMIT)));
  writeln('factor beyond LIMIT ',LIMIT_mul);

  n := 0;
  Init_Sieve(n);

  pUntouch[1] := 1;//all primes
  repeat
    k := GetNextPrimeDecomp^.pfSumOfDivs-n;
    inc(n);//n-> n+1
    if k <= LIMIT then
    begin
      If k <> 1 then
        pUntouch[k] := 1
      else
      begin
        //n-1 is prime p
        //mark p*p
        pUntouch[n] := 1;
        //mark 2*p
        //5 marked by prime 2 but that is p*p, but 4 has factor sum = 3
        pUntouch[n+2] := 1;
      end;
    end;
  until n > LIMIT-2;
  //unmark 5 and mark 0
  puntouch[5] := 0;
  pUntouch[0] := 1;

  //n=limit-1
  k := GetNextPrimeDecomp^.pfSumOfDivs-n;
  inc(n);
  If (k <> 1) AND (k<=LIMIT) then
    pUntouch[k] := 1
  else
    pUntouch[n] := 1;
  //n=limit
  k := GetNextPrimeDecomp^.pfSumOfDivs-n;
  If Not(odd(k)) AND (k<=LIMIT) then
    pUntouch[k] := 1;


  n:= limit+1;
  writeln('runtime for n<= LIMIT ',(GetTickCount64-T0)/1000:0:3,' s');
  writeln('Check the rest ',Numb2USA(IntToStr((LIMIT_mul-1)*Limit)));
  TD := GettickCount64;
  CheckRest(n,pUntouch);
  writeln('runtime ',(GetTickCount64-T0)/1000:0:3,' s');
  T0 := GetTickCount64-T0;

  OutCounts(pUntouch);
end.

Output:

TIO.RUN
LIMIT = 5,000,000
factor beyond LIMIT 171
runtime for n smaller than LIMIT 0.204 s
Check the rest 850,000,000
runtime 32.643 s
          100            5
          200           10
          300           22
          400           30
          500           38
          600           48
          700           55
          800           69
          900           80

        1,000           89
        2,000          196
        3,000          318
        4,000          443
        5,000          570
        6,000          689
        7,000          801
        8,000          936
        9,000        1,082

       10,000        1,212
       20,000        2,566
       30,000        3,923
       40,000        5,324
       50,000        6,705
       60,000        8,153
       70,000        9,586
       80,000       10,994
       90,000       12,429

      100,000       13,863
      200,000       28,572
      300,000       43,515
      400,000       58,459
      500,000       73,565
      600,000       88,828
      700,000      104,062
      800,000      119,302
      900,000      134,758

    1,000,000      150,232
    2,000,000      305,290
    3,000,000      462,110
    4,000,000      619,638
    5,000,000      777,672

Real time: 32.827 s CPU share: 99.30 %
//url=https://math.dartmouth.edu/~carlp/uupaper3.pdf
   100000    13863
   200000    28572
   300000    43515
   400000    58459
   500000    73565
   600000    88828
   700000   104062
   800000   119302
   900000   134758
  1000000   150232
  2000000   305290
  3000000   462110
  4000000   619638
  5000000   777672
  6000000   936244
  7000000  1095710
  8000000  1255016
  9000000  1414783
 10000000  1574973
 20000000  3184111
 30000000  4804331
 40000000  6430224
 50000000  8060163
 60000000  9694467
 70000000 11330312
 80000000 12967239
 90000000 14606549
100000000 16246940

... at home up to 1E8
   50,000,000    8,060,163
   60,000,000    9,694,467
   70,000,000   11,330,312
   80,000,000   12,967,239
   90,000,000   14,606,549

  100,000,000   16,246,940

real    18m51,234s

Perl

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

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 platform()!=JS and time()>t1 then
            progress("j:%,d/%,d (%3.2f%%)\r",{j,lim,(j/lim)*100})
            t1 = time()+1
        end if
    end for
    if platform()!=JS then progress("") end if
    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
    string t = elapsed(time()-t0,1," (%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-(platform()==JS))

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

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

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 j=@.#+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

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.

Version 1

Run time about 70 seconds on my Core i7 machine.

import "./math" for Int, Nums
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:")
Fmt.tprint("$,6d", untouchable.where { |n| n <= 2000 }, 10)
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

Version 2

Translation of: C

This version uses a much more efficient algorithm for calculating the sums of divisors in bulk.

Run time for untouchable numbers up to 100,000 (m = 14) is now only 1.4 seconds and 1,000,000 (m = 63) is reached in 132 seconds.

import "./math" for Int, Nums
import "./fmt" for Fmt

var limit = 1e6
var m = 63
var c = Int.primeSieve(limit, false)
var n = m * limit + 1
var sumDivs = List.filled(n, 0)
for (i in 1...n) {
    var j = i
    while (j < n) {
        sumDivs[j] = sumDivs[j] + i
        j = j + i
    }
}
var s = List.filled(n, false)
for (i in 1...n) {
    var sum = sumDivs[i] - i // proper divs sum
    if (sum <= n) s[sum] = true
}
var untouchable = [2, 5]
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:")
Fmt.tprint("$,6d", untouchable.where { |n| n <= 2000 }, 10)
System.print()
Fmt.print("$,7d 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("$,7d untouchable numbers were found  <= $,9d", count-1, p)
        p = p * 10
        if (p == limit) break
    }
}
Fmt.print("$,7d untouchable numbers were found  <= $,d", untouchable.count, limit)

Output:

As Version 1, plus:

150,232 untouchable numbers were found  <= 1,000,000