A Duffinian number is a composite number k that is relatively prime to its sigma sum σ.

The sigma sum of k is the sum of the divisors of k.



E.G.
161 is a Duffinian number.

 It is composite. (7 × 23)
 The sigma sum 192 (1 + 7 + 23 + 161) is relatively prime to 161. 


Duffinian numbers are very common.

It is not uncommon for two consecutive integers to be Duffinian (a Duffinian twin) (8, 9), (35, 36), (49, 50), etc.

Less common are Duffinian triplets; three consecutive Duffinian numbers. (63, 64, 65), (323, 324, 325), etc.

Much, much less common are Duffinian quadruplets and quintuplets. The first Duffinian quintuplet is (202605639573839041, 202605639573839042, 202605639573839043, 202605639573839044, 202605639573839045).

It is not possible to have six consecutive Duffinian numbers



Task
 Find and show the first 50 Duffinian numbers.
 Find and show at least the first 15 Duffinian triplets.


See also
 Numbers Aplenty - Duffinian numbers
 OEIS:A003624 - Duffinian numbers: composite numbers k relatively prime to sigma(k)




func IsPrime(N);        \Return 'true' if N is prime
int  N, I;
[if N <= 2 then return N = 2;
if (N&1) = 0 then \even >2\ return false;
for I:= 3 to sqrt(N) do
    [if rem(N/I) = 0 then return false;
    I:= I+1;
    ];
return true;
];
 
func    SumDiv(Num);    \Return sum of proper divisors of Num
int     Num, Div, Sum, Quot;
[Div:= 2;
Sum:= 0;
loop    [Quot:= Num/Div;
        if Div > Quot then quit;
        if rem(0) = 0 then
            [Sum:= Sum + Div;
            if Div # Quot then Sum:= Sum + Quot;
            ];
        Div:= Div+1;
        ];
return Sum+1;
];
 
func GCD(A, B);         \Return greatest common divisor of A and B
int  A, B;
[while A#B do
  if A>B then A:= A-B
         else B:= B-A;
return A;
];
 
func Duff(N);           \Return 'true' if N is a Duffinian number
int  N;
[if IsPrime(N) then return false;
return GCD(SumDiv(N), N) = 1;
];
 
int C, N;
[Format(4, 0);
C:= 0;  N:= 4;
loop    [if Duff(N) then
            [RlOut(0, float(N));
            C:= C+1;
            if C >= 50 then quit;
            if rem(C/20) = 0 then CrLf(0);
            ];
        N:= N+1;
        ];
CrLf(0);  CrLf(0);
Format(5, 0);
C:= 0;  N:= 4;
loop    [if Duff(N) & Duff(N+1) & Duff(N+2) then
            [RlOut(0, float(N));  RlOut(0, float(N+1));  RlOut(0, float(N+2));
            CrLf(0);
            C:= C+1;
            if C >= 15 then quit;
            ];
        N:= N+1;
        ];
]