56 lines
1.8 KiB
Text
56 lines
1.8 KiB
Text
divisor_sum(N) = R =>
|
|
Total = 1,
|
|
Power = 2,
|
|
% Deal with powers of 2 first
|
|
while (N mod 2 == 0)
|
|
Total := Total + Power,
|
|
Power := Power*2,
|
|
N := N div 2
|
|
end,
|
|
% Odd prime factors up to the square root
|
|
P = 3,
|
|
while (P*P =< N)
|
|
Sum = 1,
|
|
Power1 = P,
|
|
while (N mod P == 0)
|
|
Sum := Sum + Power1,
|
|
Power1 := Power1*P,
|
|
N := N div P
|
|
end,
|
|
Total := Total * Sum,
|
|
P := P+2
|
|
end,
|
|
% If n > 1 then it's prime
|
|
if N > 1 then
|
|
Total := Total*(N + 1)
|
|
end,
|
|
R = Total.
|
|
|
|
% See https://en.wikipedia.org/wiki/Aliquot_sequence
|
|
aliquot_sequence(N,Limit,Seq,Class) =>
|
|
aliquot_sequence(N,Limit,[N],Seq,Class).
|
|
|
|
aliquot_sequence(_,0,_,Seq,Class) => Seq = [], Class = 'non-terminating'.
|
|
aliquot_sequence(_,_,[0|_],Seq,Class) => Seq = [0], Class = terminating.
|
|
aliquot_sequence(N,_,[N,N|_],Seq,Class) => Seq = [], Class = perfect.
|
|
aliquot_sequence(N,_,[N,_,N|_],Seq,Class) => Seq = [N], Class = amicable.
|
|
aliquot_sequence(N,_,[N|S],Seq,Class), membchk(N,S) =>
|
|
Seq = [N], Class = sociable.
|
|
aliquot_sequence(_,_,[Term,Term|_],Seq,Class) => Seq = [], Class = aspiring.
|
|
aliquot_sequence(_,_,[Term|S],Seq,Class), membchk(Term,S) =>
|
|
Seq = [Term], Class = cyclic.
|
|
aliquot_sequence(N,Limit,[Term|S],Seq,Class) =>
|
|
Seq = [Term|Rest],
|
|
Sum = divisor_sum(Term),
|
|
Term1 is Sum - Term,
|
|
aliquot_sequence(N,Limit-1,[Term1,Term|S],Rest,Class).
|
|
|
|
main =>
|
|
foreach (N in [11,12,28,496,220,1184,12496,1264460,790,909,562,1064,1488,15355717786080,153557177860800])
|
|
aliquot_sequence(N,16,Seq,Class),
|
|
printf("%w: %w, sequence: %w ", N, Class, Seq[1]),
|
|
foreach (I in 2..len(Seq), break(Seq[I] == Seq[I-1]))
|
|
printf("%w ", Seq[I])
|
|
end,
|
|
nl
|
|
end.
|