%I #34 Sep 08 2022 08:46:19
%S 3,7,13,31,127,307,1723,2801,3541,8191,19531,86143,131071,492103,
%T 524287,552793,684757,704761,735307,797161,1353733,1886503,3413257,
%U 3894703,5473261,7094233,7781311,9250723,10378063,12655807,18143341,19443691,22292563,23907211
%N Primes q of the form sigma((p + 1) / 2) where p is a prime.
%C Mersenne primes p = 2^k - 1 (A000668) are terms: sigma((p + 1) / 2) = sigma((2^k - 1 + 1) / 2) = sigma((2^(k - 1)) = 2^k - 1 = p.
%C 2801 is the smallest term of the form 6*k + 5. The next one is 39449441. Note that both of them are of the form 1 + t + t^2 + t^3 + t^4 where t is a prime number. - _Altug Alkan_, Oct 03 2017
%e Prime 13 is a term because there is prime 17 with sigma((17 + 1) / 2) = sigma(9) = 13.
%t max = 10^6; Select[Union@ Reap[Do[If[PrimeQ@ #, Sow@ #] &@DivisorSigma[1, (Prime@ i + 1)/2], {i, max}] ][[-1, 1]], # < Prime[max]/2 &] (* _Michael De Vlieger_, Sep 16 2017, corrected by _Amiram Eldar_, Oct 08 2021 *)
%o (Magma) m := 5*10^7; Set(Sort([SumOfDivisors((n+1) div 2): n in [1..2*m] | IsPrime(n) and IsPrime(SumOfDivisors((n+1) div 2)) and SumOfDivisors((n+1) div 2) le m])) // corrected by _Amiram Eldar_, Oct 08 2021
%o (PARI) lista(nn) = {my(list = List()); forprime(p=3, 2*nn, if (isprime(q=sigma((p+1)/2)), listput(list, q));); select(x->(x <= nn), vecsort(Vec(list)));} \\ _Michel Marcus_, Oct 08 2021
%Y Cf. A000203, A000668, A292447.
%K nonn
%O 1,1
%A _Jaroslav Krizek_, Sep 16 2017