A279912.sf

#!/usr/bin/ruby

# Solved: https://oeis.org/A279912
# a(n) = Sum_{d|n} d * phi(n/d)^2
# a(p^k) = p^(k-1) * ((p-1) * p^k + 1)
# a(p^k) = p^(k-1) - p^(2k-1) + p^(2k)
# 0, 1, 3, 7, 10, 21, 21, 43, 36, 57, 63, 111, 70, 157, 129, 147, 136, 273, 171, 343, 210, 301, 333, 507, 252, 505, 471, 495, 430, 813, 441, 931, 528, 777, 819, 903, 570, 1333, 1029, 1099, 756, 1641, 903, 1807, 1110, 1197, 1521, 2163, 952, 2065, 1515, 1911, 1570, 2757

# TODO:   https://oeis.org/A279911
# 0, 1, 2, 4, 6, 11, 10, 22, 22, 31, 28, 56, 36, 79, 58, 72, 86, 137, 80, 172, 112, 145, 148, 254, 146, 261, 208, 274, 230, 407, 182, 466, 342, 375, 360, 448, 322, 667, 456, 528, 444, 821, 384, 904, 592, 635, 676, 1082, 574, 1051, 692, 924, 836, 1379, 732, 1154, 912, 1153

func f(n) {     # https://oeis.org/A279911
    sum(1..n, {|i|
      de(n**i / i)
    })
}

func f2(n) {
    map(1..n, {|i|
      de(n**i / i)
    })
}

func g2(n) {        # https://oeis.org/A279912
    map(1..n, {|i|
      de(i**n / n)
    })
}

func moebius_transform(n, f={.sigma}) {
    n.divisor_sum {|d|
        moebius(d) * f(n/d)
    }
}

func f3(n) {

    return (n.divisors.sum {|d|
        #faulhaber(n, 1) - moebius_transform(d, {_}).faulhaber(1) + 1

        #moebius_transform(d, {_}).faulhaber(1)

        #faulhaber(moebius(d) * (n/d), 1) - 1
        #d.divisors.sum {|f|
        #    moebius(f) * (d/f)
        #} * euler_phi(n/d)
    })

   n.factor_prod {|p,e|
        #(p**e * (p**e + 1) / 2 - p**e + 1)

        #(((p-1)**e * p**(e-1)).faulhaber(1)) * p**(e-1) + 1

        #faulhaber((p-1)**e, 1) + 1     # ok

        (p-1)**e * ((p-1)**e + p**(e-1)) / 2 + 1

        #p**(e-1) - p**(2*e - 1) + p**(2*e)

        #dedekind_psi((p-1).faulhaber(1) - p**e - 1)
        #dedekind_psi(p + 1)
    }
}

say 20.of(f)
say 20.of(f3)

say ''

say f(29)
say f3(29)

say ''

say f(3*5)
say f3(3*5)

__END__

say ''

say f2(3*5)
say g2(3*5)
say @(1..15)

__END__

say f2(23)
say f2(2)
say f2(3)
say f2(4)
say f2(5)
say f2(6)

say f(29)


__END__

func moebius_transform(n, f={.sigma}) {
    n.divisor_sum {|d|
        moebius(d) * f(n/d)
    }
}

func f2(n) {
    n.divisors.sum{|d|
        moebius_transform(d, { _ }).faulhaber(1)
    }
}

func g(n) {
    sum(1..n, {|i|
      de(i**n / n)
    })
}

say 20.of {|n|
    n.divisors.sum{|d|
        d * euler_phi(n/d)**2
    }
}



func inverse_moebius_transform(n, f={.sigma}) {
    n.divisor_sum {|d|
        f(d)
    }
}

say 20.of(g)

say 20.of {|n|
    n.divisors.sum{|d|
        (n/d) * moebius_transform(d, { _ })**2
    }
}

say ''

say 20.of(f)
#say {f2(_) != f(_) }.grep(1..30)

say 20.of {|n|
    n.divisors.sum{|d|
        (n/d) * moebius_transform(d, {  }).faulhaber(1)
    } - n + 1
}


__END__

say prod(1..1000, {

    ChebyshevU(3, _) / ChebyshevT(3, _)
})

say 20.of {

    ChebyshevT(3, _)
}

__END__
func f(n) {
    n.factor_prod {|p,e|
        (p/(p == 2 ? 1 : p.prev_prime))**e
    }
}

say f(100!)
say 100!.factor_exp
say f(100!).factor_exp

#say 20.of { f(_!) }

__END__

func f(n) {
    n.factor_prod{|p,k|
         #p**(k-1) * ((p+1) * p**k )
         p**(2*k + 1) * ((p - 1) )
    }
}

say 20.of(f)

__END__

func a((1)) { 1 }
func a((2)) { 3 }
func a(n) is cached { prime(n) * a(n-1) + a(n-2) }

for k in (1..5000) {
    #say [k, a(k).factor]
    if (a(k).is_prime) {
        print(k, ", ")
    }
}

#say 20.of { a(_+1) }
#say a(100)

#say 100.prime.primes.flip.reduce {|a,b| b + 1/a }.de

__END__

func E(n, m) {
    sum(1..n, {|k|
        dedekind_psi(k, m)
    })
}

func D(n, m) {
    sum(1..n, {|k|
        E(floor(n/k), m)
    })
}

func D2(n, m) {
    sum(1..n, {|k|
        k.divisors.sum {|d|
            dedekind_psi(d, m)
        }
    })
}

say 20.of {D(_, 0) }
say 20.of {D(_, 1) }
say 20.of {D(_, 2) }

say ''
say 20.of {D2(_, 0) }
say 20.of {D2(_, 1) }
say 20.of {D2(_, 2) }