login

Year-end appeal: Please make a donation to the OEIS Foundation to support ongoing development and maintenance of the OEIS. We are now in our 61st year, we have over 378,000 sequences, and we’ve reached 11,000 citations (which often say “discovered thanks to the OEIS”).

A248909
Completely multiplicative with a(p) = p if p = 6k+1 and a(p) = 1 otherwise.
8
1, 1, 1, 1, 1, 1, 7, 1, 1, 1, 1, 1, 13, 7, 1, 1, 1, 1, 19, 1, 7, 1, 1, 1, 1, 13, 1, 7, 1, 1, 31, 1, 1, 1, 7, 1, 37, 19, 13, 1, 1, 7, 43, 1, 1, 1, 1, 1, 49, 1, 1, 13, 1, 1, 1, 7, 19, 1, 1, 1, 61, 31, 7, 1, 13, 1, 67, 1, 1, 7, 1, 1, 73, 37, 1, 19, 7, 13, 79, 1, 1
OFFSET
1,7
COMMENTS
To compute a(n) replace primes not of the form 6k+1 in the prime factorization of n by 1.
The first place this sequence differs from A170824 is at n = 49.
For p prime, a(p) = p if p is a term in A002476 and a(p) = 1 if p = 2, p = 3 or p is a term in A007528.
a(n) is the largest term of A004611 that divides n. - Peter Munn, Mar 06 2021
FORMULA
a(1) = 1; for n > 1, if A020639(n) = 1 (mod 6), a(n) = A020639(n) * a(A032742(n)), otherwise a(n) = a(A028234(n)). - Antti Karttunen, Jul 09 2017
a(n) = a(A065330(n)). - Peter Munn, Mar 06 2021
EXAMPLE
a(49) = 49 because 49 = 7^2 and 7 = 6*1 + 1.
a(15) = 1 because 15 = 3*5 and neither of these primes is of the form 6k+1.
a(62) = 31 because 62 = 31*2, 31 = 6*5 + 1, and 2 is not of the form 6k+1.
MAPLE
A248909 := proc(n)
local a, pf;
a := 1 ;
for pf in ifactors(n)[2] do
if modp(op(1, pf), 6) = 1 then
a := a*op(1, pf)^op(2, pf) ;
end if;
end do:
a ;
end proc: # R. J. Mathar, Mar 14 2015
MATHEMATICA
f[p_, e_] := If[Mod[p, 6] == 1, p^e, 1]; a[1] = 1; a[n_] := Times @@ f @@@ FactorInteger[n]; Array[a, 100] (* Amiram Eldar, Sep 19 2020 *)
PROG
(Sage)
n=100; sixnplus1Primes=[x for x in primes_first_n(100) if (x-1)%6==0]
[prod([(x^(x in sixnplus1Primes))^y for x, y in factor(n)]) for n in [1..n]]
(PARI) a(n) = {my(f = factor(n)); for (i=1, #f~, if ((f[i, 1] - 1) % 6, f[i, 1] = 1); ); factorback(f); } \\ Michel Marcus, Mar 11 2015
(Python)
from sympy import factorint
def A248909(n):
y = 1
for p, e in factorint(n).items():
y *= (1 if (p-1) % 6 else p)**e
return y # Chai Wah Wu, Mar 15 2015
(Scheme) (define (A248909 n) (if (= 1 n) n (* (if (= 1 (modulo (A020639 n) 6)) (A020639 n) 1) (A248909 (A032742 n))))) ;; Antti Karttunen, Jul 09 2017
CROSSREFS
Sequences used in a definition of this sequence: A002476, A004611, A007528, A020639, A028234, A032742.
Equivalent sequence for distinct prime factors: A170824.
Equivalent sequences for prime factors of other forms: A000265 (2k+1), A343430 (3k-1), A170818 (4k+1), A097706 (4k-1), A343431 (6k-1), A065330 (6k+/-1), A065331 (<= 3).
Sequence in context: A325470 A240831 A170824 * A140213 A331927 A285483
KEYWORD
nonn,mult,easy
AUTHOR
Tom Edgar, Mar 06 2015
STATUS
approved