|
%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
|
