A267487 - OEIS

%I #13 Apr 25 2017 10:02:32

%S 2,3,7,31,127,1093,3511,8191,131071,524287

%N Primes p such that A001221(p+1)^(p-1) == 1 (mod p^2).

%C No further terms up to 10^9.

%C Are all terms of A000668 and A001220 in the sequence?

%C Does the sequence contain any terms not in A000668 or A001220 other than 2?

%p isA267487 := proc(p)

%p     if isprime(p) then

%p         A001221(p+1) ;

%p         simplify(modp(% &^ (p-1),p^2) =1 );

%p     else

%p         false;

%p     end if;

%p end proc:

%p p := 2;

%p for i from 1 do

%p     if isA267487(p) then

%p         printf("%d\n",p) ;

%p     end if;

%p     p := nextprime(p) ;

%p end do: # _R. J. Mathar_, Jan 23 2016

%t Select[Prime[Range[3200]], Mod[PrimeNu[# + 1], #^2]^(# - 1) == 1 &] (* _G. C. Greubel_, Apr 25 2017 *)

%o (PARI) forprime(p=1, 1e9, if(Mod(omega(p+1), p^2)^(p-1)==1, print1(p, ", ")))

%Y Cf. A000668, A001220, A260377.

%K nonn,more

%O 1,1

%A _Felix Fröhlich_, Jan 15 2016