login
A359062
Nonprime terms of A359059.
0
1, 8, 9, 18, 20, 27, 32, 36, 42, 44, 45, 49, 50, 54, 63, 68, 72, 78, 80, 81, 84, 90, 92, 99, 105, 108, 110, 114, 116, 117, 125, 126, 128, 135, 144, 153, 156, 162, 164, 168, 169, 170, 171, 176, 180, 186, 188, 189, 195, 198, 200, 207, 210, 212, 216, 222, 225, 228, 230
OFFSET
1,2
COMMENTS
Subsequence of A359059 after we eliminate primes.
EXAMPLE
8 is a term because 3|(4+2+12).
9 is a term because 3|(6+3+12).
18 is a term because 3|(6+8+36).
20 is a term because 3|(8+10+36).
27 is a term because 3|(18+3+36).
MAPLE
filter:= proc(n) local F, p, ph, r, ps;
F:= numtheory:-factorset(n);
if F = {n} then return false fi;
ph:= n * mul((p-1)/p, p = F);
r:= convert(F, `*`);
ps:= n * mul((p+1)/p, p = F);
(ph+r+ps) mod 3 = 0
end proc:
select(filter, [$1..1000]); # Robert Israel, Dec 20 2022
MATHEMATICA
q[n_] := Module[{f = FactorInteger[n], p, e}, p = f[[;; , 1]]; e = f[[;; , 2]]; Divisible[Times @@ ((p - 1)*p^(e - 1)) + Times @@ p + Times @@ ((p + 1)*p^(e - 1)), 3]]; Select[Range[230], ! PrimeQ[#] && q[#] &] (* Amiram Eldar, Dec 20 2022 *)
PROG
(Python)
from sympy.ntheory.factor_ import totient
from sympy import isprime, primefactors, prod
def rad(n): return 1 if n < 2 else prod(primefactors(n))
def psi(n):
plist = primefactors(n)
return n*prod(p+1 for p in plist)//prod(plist)
# Output display terms.
for n in range(1, 231):
if(False == isprime(n)):
if(0 == (totient(n) + rad(n) + psi(n)) % 3):
print(n, end = ", ")
(PARI) isok(m) = !isprime(m) && (((eulerphi(m) + factorback(factorint(m)[, 1]) + m*sumdiv(m, d, moebius(d)^2/d)) % 3) == 0); \\ Michel Marcus, Dec 27 2022
CROSSREFS
Cf. A000010 (phi), A001615 (psi), A007947 (rad).
Subsequence of A359059.
Sequence in context: A041551 A356699 A057104 * A095191 A050706 A159836
KEYWORD
nonn
AUTHOR
Torlach Rush, Dec 15 2022
STATUS
approved