%I #20 Apr 20 2024 09:57:29
%S 21,57,133,381,553,813,993,1057,1333,1561,1641,1893,1981,2653,2757,
%T 3193,3661,5257,5853,6973,8373,8557,9121,9313,10713,10921,12657,13341,
%U 15253,15501,16257,18633,19741,22053,24493,29413,30801,32221,32581,33673,35157,39801
%N Intersection of A002061 and A016105.
%C If p is a cuban prime (A002407) and p == 3 (mod 4) (A002145), then m = 3*p is a term. Indeed, there is k for which p = 1 + 3*k*(k + 1) and m = 3*p = 3 + 9*k*(k + 1) = (3*k + 2)^2 - (3*k + 2) + 1, so m is a term.
%C The sequence also includes terms that do not have this form: 133 = 12^2 - 12 + 1 = 7*19, 553 = 24^2 - 24 + 1 = 7*79, 1057 = 33^2 - 33 + 1 = 7*151, 1333 = 37^2 - 37 + 1= 31*43 and others.
%e A002061(5) = 21 = A016105(1), so 21 is a term.
%e A002061(8) = 57 = A016105(3), so 57 is a term.
%t TR=40000; R1=Ceiling[(1+Sqrt[1-4(1-TR)])/2]; R2=TR/4; Intersection[Table[n^2-n+1, {n, 0, R1}], Select[4Range[5, R2]+1, PrimeNu[#]==2&&MoebiusMu[#]==1&&Mod[FactorInteger[#][[1, 1]], 4]!=1&]](* _James C. McMahon_, Feb 27 2024 *)
%o (Magma) pd:=PrimeDivisors; blum:=func<n|#Divisors(n) eq 4 and #pd(n) eq 2 and pd(n)[1] mod 4 eq 3 and pd(n)[2] mod 4 eq 3>; [n:n in [s^2-s+1:s in [2..2000]]|blum(n)];
%Y Cf. A002061, A002145, A002407, A016105.
%K nonn
%O 1,1
%A _Marius A. Burtea_, Feb 27 2024
|