The OEIS is supported by the many generous donors to the OEIS Foundation.
%I #34 Jan 30 2025 17:22:07
%S 2,3,5,7,8,11,12,13,17,18,19,20,23,27,28,29,31,32,37,41,43,44,45,47,
%T 48,50,52,53,59,61,63,67,68,71,72,73,75,76,79,80,83,89,92,97,98,99,
%U 101,103,107,108,109,112,113,116,117,124,125,127,128,131,137,139,147,148,149
%N Numbers of the form p * m^2, where p is prime and m > 0: union of A228056 and A000040.
%C No term is the product of two other terms.
%C Squares of terms and pairwise products of distinct terms form a subsequence of A028260.
%C Numbers n such that A162642(n) = 1. - _Jason Kimberley_, Oct 10 2016
%C Numbers k such that A007913(k) is a prime number. - _Amiram Eldar_, Jul 27 2020
%H Chris Boyd, <a href="/A229125/b229125.txt">Table of n, a(n) for n = 1..10000</a>
%H Eckford Cohen, <a href="https://doi.org/10.4064/aa-7-4-417-420">Arithmetical notes, IX. On the set of integers representable as a product of a prime and square</a>, Acta Arithmetica, Vol. 7 (1962), pp. 417-420.
%F The number of terms not exceeding x is (Pi^2/6) * x/log(x) + O(x/(log(x))^2) (Cohen, 1962). - _Amiram Eldar_, Jul 27 2020
%t With[{nn=70},Take[Union[Flatten[Table[p*m^2,{p,Prime[Range[nn]]},{m,nn}]]], nn]] (* _Harvey P. Dale_, Dec 02 2014 *)
%o (PARI) test(n)=isprime(core(n))
%o for(n=1,200,if(test(n), print1(n",")))
%o (Python)
%o from math import isqrt
%o from sympy import primepi
%o def A229125(n):
%o def bisection(f,kmin=0,kmax=1):
%o while f(kmax) > kmax: kmax <<= 1
%o kmin = kmax >> 1
%o while kmax-kmin > 1:
%o kmid = kmax+kmin>>1
%o if f(kmid) <= kmid:
%o kmax = kmid
%o else:
%o kmin = kmid
%o return kmax
%o def f(x): return n+x-sum(primepi(x//y**2) for y in range(1,isqrt(x)+1))
%o return bisection(f,n,n) # _Chai Wah Wu_, Jan 30 2025
%Y Subsequence of A026424.
%Y Cf. A028260, A162642, A229153.
%Y Cf. A007913, A013661.
%Y Subsequences: A000040, A030078, A050997, A054753, A092759, A179643, A179665, A246551.
%K nonn,changed
%O 1,1
%A _Chris Boyd_, Sep 14 2013