%I #77 Feb 01 2025 12:03:09
%S 1,2,3,5,6,7,9,10,11,13,14,15,17,18,19,21,22,23,25,26,29,30,31,33,34,
%T 35,37,38,39,41,42,43,45,46,47,49,50,51,53,55,57,58,59,61,62,63,65,66,
%U 67,69,70,71,73,74,75,77,78,79,82,83,85,86,87,89,90,91,93,94,95,97,98
%N Numbers not divisible by p^p for any prime p.
%C If a(n) = Product p_i^e_i then p_i > e_i for all i.
%C Complement of A100716; A129251(a(n)) = 0. - _Reinhard Zumkeller_, Apr 07 2007
%C Density is 0.72199023441955... = Product_{p>=2} (1 - p^-p) where p runs over the primes. - _Charles R Greathouse IV_, Jan 25 2012
%C A027748(a(n),k) <= A124010(a(n),k), 1<=k<=A001221(a(n)). - _Reinhard Zumkeller_, Apr 28 2012
%C Range of A276086. Also numbers not divisible by m^m for any natural number m > 1. - _Antti Karttunen_, Nov 18 2024
%H Reinhard Zumkeller, <a href="/A048103/b048103.txt">Table of n, a(n) for n = 1..10000</a>
%F a(n) ~ kn with k = 1/Product_{p>=2}(1 - p^-p) = Product_{p>=2}(1 + 1/(p^p - 1)) = 1.3850602852..., where the product is over all primes p. - _Charles R Greathouse IV_, Jan 25 2012
%F For n >= 1, A377982(a(n)) = n. - _Antti Karttunen_, Nov 18 2024
%e 6 = 2^1 * 3^1 is OK but 12 = 2^2 * 3^1 is not.
%e 625 = 5^4 is present because it is not divisible by 5^5.
%t {1}~Join~Select[Range@ 120, Times @@ Boole@ Map[First@ # > Last@ # &, FactorInteger@ #] > 0 &] (* _Michael De Vlieger_, Aug 19 2016 *)
%o (Haskell)
%o a048103 n = a048103_list !! (n-1)
%o a048103_list = filter (\x -> and $
%o zipWith (>) (a027748_row x) (map toInteger $ a124010_row x)) [1..]
%o -- _Reinhard Zumkeller_, Apr 28 2012
%o (Scheme, with Antti Karttunen's IntSeq-library)
%o (define A048103 (ZERO-POS 1 1 A129251))
%o ;; _Antti Karttunen_, Aug 18 2016
%o (PARI) isok(n) = my(f=factor(n)); for (i=1, #f~, if (f[i,1] <= f[i,2], return(0))); return(1); \\ _Michel Marcus_, Nov 13 2020
%o (PARI) A359550(n) = { my(pp); forprime(p=2, , pp = p^p; if(!(n%pp), return(0)); if(pp > n, return(1))); }; \\ (A359550 is the characteristic function for A048103) - _Antti Karttunen_, Nov 18 2024
%o (Python)
%o from itertools import count, islice
%o from sympy import factorint
%o def A048103_gen(startvalue=1): # generator of terms >= startvalue
%o return filter(lambda n:all(map(lambda d:d[1]<d[0],factorint(n).items())),count(max(startvalue,1)))
%o A048103_list = list(islice(A048103_gen(),30)) # _Chai Wah Wu_, Jan 05 2023
%Y Complement: A100716.
%Y Positions of 0's in A129251, A342023, A376418, positions of 1's in A327936, A342007, A359550 (characteristic function).
%Y Cf. A048102, A048104, A051674 (p^p), A054743, A054744, A377982 (a left inverse, partial sums of char. fun, see also A328402).
%Y Cf. A276086 (permutation of this sequence, see also A376411, A376413).
%Y Subsequences: A002110, A005117, A006862, A024451 (after its initial 0), A057588, A099308 (after its initial 0), A276092, A328832, A359547, A370114, A371083, A373848, A377871, A377992.
%Y Disjoint union of {1}, A327934 and A358215.
%Y Also A276078 is a subsequence, from which this differs for the first time at n=451 where a(451)=625, while that value is missing from A276078.
%K nonn,easy
%O 1,2
%A _N. J. A. Sloane_
%E More terms from _James A. Sellers_, Apr 22 2000