This is a front-end for the Online Encyclopedia of Integer Sequences, made by Christian Perfect. The idea is to provide OEIS entries in non-ancient HTML, and then to think about how they're presented visually. The source code is on GitHub.
%I A370519 #27 Feb 28 2025 12:06:51 %S A370519 21,57,133,381,553,813,993,1057,1333,1561,1641,1893,1981,2653,2757, %T A370519 3193,3661,5257,5853,6973,8373,8557,9121,9313,10713,10921,12657,13341, %U A370519 15253,15501,16257,18633,19741,22053,24493,29413,30801,32221,32581,33673,35157,39801 %N A370519 Intersection of A002061 and A016105. %C A370519 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 A370519 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. %H A370519 Robert Israel, <a href="/A370519/b370519.txt">Table of n, a(n) for n = 1..572</a> %e A370519 A002061(5) = 21 = A016105(1), so 21 is a term. %e A370519 A002061(8) = 57 = A016105(3), so 57 is a term. %p A370519 N:= 10^5: # for terms <= N %p A370519 P:= select(isprime, [seq(i,i=3..N/3,4)]): %p A370519 sort(select(t -> t <= N and issqr(4*t-3), [seq(seq(P[i]*P[j],i=1..j-1),j=1..nops(P))])); # _Robert Israel_, Feb 27 2025 %t A370519 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 A370519 (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 A370519 Cf. A002061, A002145, A002407, A016105. %K A370519 nonn %O A370519 1,1 %A A370519 _Marius A. Burtea_, Feb 27 2024