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.