login
a(n) is the product of primes (with multiplicity) of form 4*k+1 that divide n.
10

%I #41 Jan 05 2025 19:51:39

%S 1,1,1,1,5,1,1,1,1,5,1,1,13,1,5,1,17,1,1,5,1,1,1,1,25,13,1,1,29,5,1,1,

%T 1,17,5,1,37,1,13,5,41,1,1,1,5,1,1,1,1,25,17,13,53,1,5,1,1,29,1,5,61,

%U 1,1,1,65,1,1,17,1,5,1,1,73,37,25,1,1,13,1,5,1,41,1,1,85,1,29,1

%N a(n) is the product of primes (with multiplicity) of form 4*k+1 that divide n.

%C Completely multiplicative with a(p) = p if p = 4k+1 and a(p) = 1 otherwise. - _Tom Edgar_, Mar 05 2015

%H Alois P. Heinz, <a href="/A170818/b170818.txt">Table of n, a(n) for n = 1..10000</a>

%H A. Tripathi, <a href="https://web.archive.org/web/2024*/https://www.fq.math.ca/Papers1/46_47-4/Tripathi.pdf">On Pythagorean triples containing a fixed integer</a>, Fib. Q., 46/47 (2008/2009), 331-340.

%H <a href="/index/Di#divseq">Index to divisibility sequences</a>

%F a(n) = n/A072438(n). - _Michel Marcus_, Mar 05 2015

%p a:= n-> mul(`if`(irem(i[1], 4)=1, i[1]^i[2], 1), i=ifactors(n)[2]):

%p seq(a(n), n=1..100); # _Alois P. Heinz_, Jun 09 2014

%t a[n_] := Product[{p, e} = pe; If[Mod[p, 4] == 1, p^e, 1], {pe, FactorInteger[n]}];

%t Array[a, 100] (* _Jean-François Alcover_, May 29 2019 *)

%o (PARI) a(n)=my(f=factor(n)); prod(i=1,#f~,if(f[i,1]%4>1,1,f[i,1])^f[i,2]) \\ _Charles R Greathouse IV_, Jun 28 2015

%o (Python)

%o from sympy import factorint, prod

%o def a072438(n):

%o f = factorint(n)

%o return 1 if n == 1 else prod(i**f[i] for i in f if i % 4 != 1)

%o def a(n): return n//a072438(n) # _Indranil Ghosh_, May 08 2017

%Y Cf. A170817-A170819, A097706, A083025, A072438, A286361.

%K nonn,mult

%O 1,5

%A _N. J. A. Sloane_, Dec 22 2009