A391004 Primes p that can be represented as the sum of distinct nontrivial divisors of p + 1.

5, 11, 17, 19, 23, 29, 41, 47, 53, 59, 71, 79, 83, 89, 103, 107, 131, 139, 149, 167, 179, 191, 197, 199, 223, 227, 233, 239, 251, 263, 269, 271, 293, 307, 311, 349, 359, 367, 379, 383, 389, 419, 431, 439, 449, 461, 463, 467, 479, 499, 503, 509, 557, 569, 571, 587

Offset

1,1

Comments

Nontrivial divisors of n are > 1 and < n (A070824).

Subsequence of A391005 (semiperfect primes).

In many cases p + 1 is a practical number (A005153). In fact in the above list only 349 + 1 and 571 + 1 are not practical numbers.

Links

Peter Luschny, Table of n, a(n) for n = 1..10000

Peter Luschny, Divisor subsets. A Python implementation for A39100[4/6/7].

Example

   p as sum of subset of divisors of p + 1.
  --------------------------------------------
   5 = 2 +  3;
  11 = 2 +  3 +  6;
  17 = 2 +  6 +  9;
  19 = 4 +  5 + 10;
  23 = 3 +  8 + 12;
  29 = 3 +  5 +  6 + 15;
  41 = 6 + 14 + 21;
  47 = 3 +  8 + 12 + 24;
  53 = 2 +  6 + 18 + 27;
  59 = 4 + 10 + 15 + 30;
  71 = 8 +  9 + 18 + 36;
  79 = 5 +  8 + 10 + 16 + 40;
  499= 4 + 20 + 100 + 125 + 250;
.
347 is not a term, even though it is referenced in A391003, because A391003 does not require the divisors to be nontrivial.

Maple

issum:= proc(x, S)
  option remember;
  local Sp, s;
  Sp:= select(`<=`, S, x);
  if member(x, Sp) then return true fi;
  if convert(Sp, `+`) < x then return false fi;
  s:= max(Sp); Sp:= Sp minus {s};
  procname(x-s, Sp) or procname(x, Sp)
end proc:  # Robert Israel, Dec 01 2025
filter:= proc(n)
  local S;
  S:= numtheory:-divisors(n+1) minus {1, n+1};
  issum(n, S)
end proc:
select(filter, {seq(ithprime(i), i=1..1229)});

Mathematica

Select[Prime[Range[3, 100]], MemberQ[Total/@Subsets[Rest[Drop[Divisors[#+1], -1]]], #]&] (* James C. McMahon, Dec 01 2025 *)

Prog

(SageMath)
def is_A391004(n):
    if not is_prime(n): return False
    divs = sorted(divisors(n + 1))[:-1]
    k = len(divs)
    for mask in range(1, 1 << k):
        subset = [divs[i] for i in range(k) if mask & (1 << i)]
        if sum(subset) == n + 1 and 1 in subset: return True
    return False
print([n for n in range(2, 600) if is_A391004(n)])
(PARI) is_a391004(n) = my(m=n+1, s=0, so); if(!isprime(n), return(0)); my(dnt=setminus(divisors(m), [1, m]), dodd=select(x->x%2, dnt), nodd=#dodd, deven=setminus(dnt, dodd), neven=#deven); if(nodd==0, return(0)); forstep(jodd=1, nodd, 2, forsubset([nodd, jodd], sodd, so=sum(k=1, #sodd, dodd[sodd[k]]); if(so==n, return(1)); forsubset(neven, seven, my(ns=#seven); if(ns>0, s = so + sum(k=1, ns, deven[seven[k]])); if(s==n, return(1))))); 0 \\ Hugo Pfoertner, Nov 29 2025

Crossrefs

Cf. A000040, A070824, A391005, A391003, A070824, A085498, A085499, A005153.

Sequence in context: A096264 A391005 A391003 * A239709 A175249 A271982

Adjacent sequences: A391001 A391002 A391003 * A391005 A391006 A391007

Keyword

nonn

Author

Peter Luschny, Nov 25 2025

Status

approved