%I #20 Mar 08 2024 08:08:51
%S 1,5,13,17,29,37,41,53,61,65,73,85,89,97,101,109,113,137,145,149,157,
%T 173,181,185,193,197,205,221,229,233,241,257,265,269,277,281,293,305,
%U 313,317,337,349,353,365,373,377,389,397,401,409,421,433,445,449
%N Products of distinct primes congruent to 1 modulo 4 (A002144).
%C Contains A002144 as a subsequence, and is a subsequence of A016813 and of A005117.
%C Also, these numbers satisfy A231589(n) = floor(n*(n-1)/4) (A011848).
%H R. J. Mathar, <a href="/A231754/b231754.txt">Table of n, a(n) for n = 1..4999</a>
%H Rafael Jakimczuk, <a href="http://dx.doi.org/10.13140/RG.2.2.27745.48487">Generalizations of Mertens's Formula and k-Free and s-Full Numbers with Prime Divisors in Arithmetic Progression</a>, ResearchGate, 2024.
%H Shailesh A. Shirali, <a href="http://www.jstor.org/stable/2690862">A family portrait of primes-a case study in discrimination</a>, Math. Mag., Vol. 70, No. 4 (Oct., 1997), pp. 263-272.
%F The number of terms that do not exceed x is ~ c * x / sqrt(log(x)), where c = A088539 * sqrt(A175647) / Pi = 0.3097281805... (Jakimczuk, 2024, Theorem 3.10, p. 26). - _Amiram Eldar_, Mar 08 2024
%e 65 = 5*13 is in the sequence since both 5 and 13 are congruent to 1 modulo 4.
%p isA231754 := proc(n)
%p local d;
%p for d in ifactors(n)[2] do
%p if op(2,d) > 1 then
%p return false;
%p elif modp(op(1,d),4) <> 1 then
%p return false;
%p end if;
%p end do:
%p true ;
%p end proc:
%p for n from 1 to 500 do
%p if isA231754(n) then
%p printf("%d,",n) ;
%p end if;
%p end do: # _R. J. Mathar_, Mar 16 2016
%t Select[Range[500], # == 1 || AllTrue[FactorInteger[#], Last[#1] == 1 && Mod[First[#1], 4] == 1 &] &] (* _Amiram Eldar_, Mar 08 2024 *)
%o (PARI) isok(n) = if (! issquarefree(n), return (0)); if (n > 1, f = factor(n); for (i=1, #f~, if (f[i, 1] % 4 != 1, return (0)))); 1
%Y Intersection of A005117 and A004613.
%Y Cf. A002144, A011848, A016813, A167181, A231589.
%Y Cf. A088539, A175647.
%K nonn,easy
%O 1,2
%A _Michel Marcus_, Nov 13 2013