A103339 Numerator of the unitary harmonic mean (i.e., the harmonic mean of the unitary divisors) of the positive integer n.
1, 4, 3, 8, 5, 2, 7, 16, 9, 20, 11, 12, 13, 7, 5, 32, 17, 12, 19, 8, 21, 22, 23, 8, 25, 52, 27, 14, 29, 10, 31, 64, 11, 68, 35, 72, 37, 38, 39, 80, 41, 7, 43, 44, 3, 23, 47, 48, 49, 100, 17, 104, 53, 18, 55, 28, 57, 116, 59, 4, 61, 31, 63, 128, 65, 11, 67, 136, 23, 35, 71, 16, 73
Offset: 1
Examples
1, 4/3, 3/2, 8/5, 5/3, 2, ...
a(8) = 16 because the unitary divisors of 8 are {1,8} and 2/(1/1 + 1/8) = 16/9.
Links
- Reinhard Zumkeller, Table of n, a(n) for n = 1..10000
Programs
-
Haskell
import Data.Ratio ((%), numerator) a103339 = numerator . uhm where uhm n = (n * a034444 n) % (a034448 n) -- Reinhard Zumkeller, Mar 17 2012
-
Maple
with(numtheory): udivisors:=proc(n) local A, k: A:={}: for k from 1 to tau(n) do if gcd(divisors(n)[k],n/divisors(n)[k])=1 then A:=A union {divisors(n)[k]} else A:=A fi od end: utau:=n->nops(udivisors(n)): usigma:=n->sum(udivisors(n)[j],j=1..nops(udivisors(n))): uH:=n->n*utau(n)/usigma(n):seq(numer(uH(n)),n=1..81); -
Mathematica
ud[n_] := 2^PrimeNu[n]; usigma[n_] := DivisorSum[n, If[GCD[#, n/#] == 1, #, 0]&]; a[1] = 1; a[n_] := Numerator[n*ud[n]/usigma[n]]; Array[a, 100] (* Jean-François Alcover, Dec 03 2016 *) a[n_] := Numerator[n * Times @@ (2 / (1 + Power @@@ FactorInteger[n]))]; a[1] = 1; Array[a, 100] (* Amiram Eldar, Mar 10 2023 *)
-
PARI
a(n) = {my(f = factor(n)); numerator(n * prod(i=1, #f~, 2/(1 + f[i, 1]^f[i, 2]))); } \\ Amiram Eldar, Mar 10 2023 -
Python
from sympy import gcd from sympy.ntheory.factor_ import udivisor_sigma def A103339(n): return (lambda x, y: y*n//gcd(x,y*n))(udivisor_sigma(n),udivisor_sigma(n,0)) # Chai Wah Wu, Oct 20 2021