%I #69 Sep 26 2024 21:24:05
%S 21,33,57,69,77,93,129,133,141,161,177,201,209,213,217,237,249,253,
%T 301,309,321,329,341,381,393,413,417,437,453,469,473,489,497,501,517,
%U 537,553,573,581,589,597,633,649,669,681,713,717,721,737,749,753,781,789
%N Blum integers: numbers of the form p * q where p and q are distinct primes congruent to 3 (mod 4).
%C Subsequence of A084109. - _Ralf Stephan_ and _David W. Wilson_, Apr 17 2005
%C Subsequence of A046388. - _Altug Alkan_, Dec 10 2015
%C Subsequence of A339817. No common terms with A339870. - _Antti Karttunen_, Dec 26 2020
%C Named after the Venezuelan-American computer scientist Manuel Blum (b. 1938). - _Amiram Eldar_, Jun 06 2021
%C First introduced by Blum, Blum, & Shub for the generation of pseudorandom numbers and later applied (by Manuel Blum and other authors) to zero-knowledge proofs. - _Charles R Greathouse IV_, Sep 26 2024
%D Lenore Blum, Manuel Blum, and Mike Shub. A simple unpredictable pseudorandom number generator, SIAM Journal on computing 15:2 (1986), pp. 364-383.
%H Antti Karttunen, <a href="/A016105/b016105.txt">Table of n, a(n) for n = 1..26828</a> (all terms < 2^19; first 1000 terms from T. D. Noe)
%H Joe Hurd, <a href="http://www.gilith.com/research/talks/cambridge1997.pdf">Blum Integers</a>, Talk at the Trinity College, Jan 20 1997.
%H Wikipedia, <a href="http://en.wikipedia.org/wiki/Blum_integer">Blum integer</a>.
%F a(n) = A195758(n) * A195759(n). - _Reinhard Zumkeller_, Sep 23 2011
%F a(n) ~ 4n log n/log log n. - _Charles R Greathouse IV_, Sep 17 2022
%p N:= 10000: # to get all terms <= N
%p Primes:= select(isprime, [seq(i,i=3..N/3,4)]):
%p S:=select(`<=`,{seq(seq(Primes[i]*Primes[j],i=1..j-1),j=2..nops(Primes))},N):
%p sort(convert(S,list)); # _Robert Israel_, Dec 11 2015
%t With[{upto = 820}, Select[Union[Times@@@Subsets[ Select[Prime[Range[ PrimePi[ NextPrime[upto/3]]]], Mod[#, 4] == 3 &], {2}]], # <= upto &]] (* _Harvey P. Dale_, Aug 19 2011 *)
%t Select[4Range[5, 197] + 1, PrimeNu[#] == 2 && MoebiusMu[#] == 1 && Mod[FactorInteger[#][[1, 1]], 4] != 1 &] (* _Alonso del Arte_, Nov 18 2015 *)
%o (Haskell)
%o import Data.Set (singleton, fromList, deleteFindMin, union)
%o a016105 n = a016105_list !! (n-1)
%o a016105_list = f [3,7] (drop 2 a002145_list) 21 (singleton 21) where
%o f qs (p:p':ps) t s
%o | m < t = m : f qs (p:p':ps) t s'
%o | otherwise = m : f (p:qs) (p':ps) t' (s' `union` (fromList pqs))
%o where (m,s') = deleteFindMin s
%o t' = head $ dropWhile (> 3*p') pqs
%o pqs = map (p *) qs
%o -- _Reinhard Zumkeller_, Sep 23 2011
%o (Perl) use ntheory ":all"; forcomposites { say if ($_ % 4) == 1 && is_square_free($_) && scalar(factor($_)) == 2 && !scalar(grep { ($_ % 4) != 3 } factor($_)); } 10000; # _Dana Jacobsen_, Dec 10 2015
%o (PARI) list(lim)=my(P=List(),v=List(),t,p); forprimestep(p=3,lim\3,4, listput(P,p)); for(i=2,#P, p=P[i]; for(j=1,i-1, t=p*P[j]; if(t>lim, break); listput(v,t))); Set(v) \\ _Charles R Greathouse IV_, Jul 01 2016, updated Sep 26 2024
%o (PARI) isA016105(n) = (2==omega(n)&&2==bigomega(n)&&1==(n%4)&&3==((factor(n)[1,1])%4)); \\ _Antti Karttunen_, Dec 26 2020
%o (Python)
%o from sympy import factorint
%o def ok(n):
%o fn = factorint(n)
%o return len(fn) == sum(fn.values()) == 2 and all(f%4 == 3 for f in fn)
%o print([k for k in range(790) if ok(k)]) # _Michael S. Branicky_, Dec 20 2021
%Y Cf. A002145, A006881, A046388, A339870.
%Y Intersection of A005117 and A107978.
%Y Also, subsequence of the following sequences: A046388, A084109, A091113, A167181, A339817.
%K nonn,easy,nice
%O 1,1
%A _Robert G. Wilson v_
%E More terms from _Erich Friedman_