login
A246544
Consider the aliquot parts, in ascending order, of a composite number. Take their sum and repeat the process deleting the minimum number and adding the previous sum. The sequence lists the numbers that after some iterations reach a sum equal to themselves.
16
6, 21, 28, 85, 496, 2133, 8128, 19521, 77125, 97273, 176661, 615281, 4948133, 33550336, 68353213, 129127041, 8589869056
OFFSET
1,1
COMMENTS
Similar to Keith numbers and Primonacci numbers, using proper divisors instead of digits or prime factors.
The perfect numbers A000396 are a subset.
The numbers of iterations are: 1, 2, 1, 3, 1, 2, 1, 2, 3, 4, 2, 5, ...; with 1 when the term is perfect. - Michel Marcus, Aug 30 2014
The 2-hyperperfect numbers A007593 are a subset, with number of iterations 2. - Michel Marcus, Sep 22 2014
EXAMPLE
Aliquot parts of 85 are 1, 5 and 17:
1 + 5 + 17 = 23;
5 + 17 + 23 = 45;
17 + 23 + 45 = 85.
Aliquot parts of 19521 are 1, 3, 9, 27, 81, 241, 723, 2169 and 6507:
1 + 3 + 9 + 27 + 81 + 241 + 723 + 2169 + 6507 = 9761;
3 + 9 + 27 + 81 + 241 + 723 + 2169 + 6507 + 9761 = 19521.
MAPLE
with(numtheory): P:=proc(q, h)
local a, b, k, n, t, v; v:=array(1..h);
for n from 2 to q do if not isprime(n) then
a:=sort([op(divisors(n))]); b:=nops(a)-1;
for k from 1 to b do v[k]:=a[k]; od;
t:=b+1; v[t]:=add(v[k], k=1..b);
while v[t]<n do t:=t+1; v[t]:=add(v[k], k=t-b..t-1); od;
if v[t]=n then print(n); fi; fi; od; end: P(10^9, 1000);
MATHEMATICA
A246544 = {};
For[n = 4, n <= 1000000, n++,
If[PrimeQ[n], Continue[]];
a = Most[Divisors[n]];
sum = Total[a];
While[sum < n, sum = Total[a = Join[Rest[a], {sum}]]];
If[sum == n, AppendTo[A246544, n]];
]; A246544 (* Robert Price, Sep 08 2019 *)
PROG
(PARI) lista(nn) = {forcomposite(n=1, nn, d = divisors(n); v = vector(#d-1, i, d[i]); vs = sum(i=1, #v, v[i]); ind = 1; while (vs < n, v = concat(v, vs); vs += vs - v[ind]; ind++; ); if (vs == n, print1(n, ", ")); ); } \\ Michel Marcus, Aug 29 2014
(Python)
import math
.
def divs(n):
....large_divisors = []
....for i in range(1, int(math.sqrt(n) + 1)):
........if n % i is 0:
............yield i
............if i is not n / i:
................large_divisors.insert(0, n / i)
....for divisor in large_divisors:
........yield divisor
.
a = 2
while a < 1000000000:
....q = list(divs(a))[:-1]
....r = sum(q)
....if r > a or len(q) == 1:
........pass
....elif r == a:
........print(a)
....else:
........c = 1
........while r < a:
............q.append(r)
............r = sum(q[c:])
............c += 1
........if r == a:
............print(a)
....a += 1
# David Consiglio, Jr., Sep 09 2014
(Python)
from sympy import divisors, isprime
A246544_list = []
for n in range(2, 10**5):
....if not isprime(n):
........x = divisors(n)
........x.pop()
........y = sum(x)
........while y < n:
............x, y = x[1:]+[y], 2*y-x[0]
........if y == n:
............A246544_list.append(n)
# Chai Wah Wu, Nov 03 2014
CROSSREFS
KEYWORD
nonn,more
AUTHOR
Paolo P. Lava, Aug 29 2014
EXTENSIONS
a(13)-a(15) from Michel Marcus, Aug 29 2014
a(16) from David Consiglio, Jr., Sep 06 2014
a(17) from Lars Blomberg, Oct 27 2014
STATUS
approved