The Dirichlet convolution is a special kind of convolution that appears as a very useful tool in number theory. It operates on the set of arithmetic functions.

Challenge

Given two arithmetic functions $f,g$ (i.e. functions $f,g: mathbb N to mathbb R$) compute the Dirichlet convolution $(f * g): mathbb N to mathbb R$ as defined below.

Details

• We use the convention $ 0 notin mathbb N = 1,2,3,ldots $.


• The Dirichlet convolution $f*g$ of two arithmetic functions $f,g$ is again an arithmetic function, and it is defined as $$(f * g)(n) = sum_limitsd fleft(fracndright)cdot g(d) = sum_icdot j = n f(i)cdot g(j).$$ (Both sums are equivalent. The expression $d|n$ means $d in mathbb N$ divides $n$, therefore the summation is over the natural divisors of $n$. Similarly we can subsitute $ i = fracnd in mathbb N, j =d in mathbb N $ and we get the second equivalent formulation. If you're not used to this notation there is a step by step example at below.) Just to elaborate (this is not directly relevant for this challenge): The definition comes from computing the product of Dirichlet series: $$left(sum_ninmathbb Nfracf(n)n^sright)cdot left(sum_ninmathbb Nfracg(n)n^sright) = sum_ninmathbb Nfrac(f * g)(n)n^s$$
  


• The input is given as two black box functions. Alternatively, you could also use an infinite list, a generator, a stream or something similar that could produce an unlimited number of values.


• There are two output methods: Either a function $f*g$ is returned, or alternatively you can take take an additional input $n in mathbb N$ and return $(f*g)(n)$ directly.


• For simplicity you can assume that every element of $ mathbb N$ can be represented with e.g. a positive 32-bit int.


• For simplicity you can also assume that every entry $ mathbb R $ can be represented by e.g. a single real floating point number.

Examples

Let us first define a few functions. Note that the list of numbers below each definition represents the first few values of that function.

• the multiplicative identity (A000007)
  $$epsilon(n) = begincases1 & n=1 \ 0 & n>1 endcases$$
  1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, ...
  


• the constant unit function (A000012)$$ mathbb 1(n) = 1 : forall n $$
  1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, ...
  


• the identity function (A000027)
  $$ id(n) = n : forall n $$
  1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, ...
  


• the Möbius function (A008683)
  $$ mu(n) = begincases (-1)^k & text if n text is squarefree and k text is the number of Primefactors of n \ 0 & text otherwise endcases $$
  1, -1, -1, 0, -1, 1, -1, 0, 0, 1, -1, 0, -1, 1, 1, 0, -1, 0, -1, ...
  


• the Euler totient function (A000010)
  $$ varphi(n) = nprod_n left( 1 - frac1pright) $$
  1, 1, 2, 2, 4, 2, 6, 4, 6, 4, 10, 4, 12, 6, 8, 8, 16, 6, 18, 8, ...
  


• the Liouville function (A008836)
  $$ lambda (n) = (-1)^k $$ where $k$ is the number of prime factors of $n$ counted with multiplicity
  1, -1, -1, 1, -1, 1, -1, -1, 1, 1, -1, -1, -1, 1, 1, 1, -1, -1, -1, ...
  


• the divisor sum function (A000203)
  $$sigma(n) = sum_d d $$
  1, 3, 4, 7, 6, 12, 8, 15, 13, 18, 12, 28, 14, 24, 24, 31, 18, 39, 20, ...
  


• the divisor counting function (A000005)
  $$tau(n) = sum_d 1 $$
  1, 2, 2, 3, 2, 4, 2, 4, 3, 4, 2, 6, 2, 4, 4, 5, 2, 6, 2, 6, 4, 4, 2, 8, ...
  


• the characteristic function of square numbers (A010052)
  $$sq(n) = begincases 1 & text if n text is a square number \ 0 & textotherwiseendcases$$
  1, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, ...
  

Then we have following examples:

• $ epsilon = mathbb 1 * mu $


• $ f = epsilon * f : forall f $


• $ epsilon = lambda * vert mu vert $


• 
  $ sigma = varphi * tau $


• 
  $ id = sigma * mu$ and $ sigma = id * mathbb 1$
  


• 
  $ sq = lambda * mathbb 1 $ and $ lambda = mu * sq$
  


• 
  $ tau = mathbb 1 * mathbb 1$ and $ mathbb 1 = tau * mu $
  


• 
  $ id = varphi * mathbb 1 $ and $ varphi = id * mu $
  

The last for are a consequence of the Möbius inversion: For any $f,g$ the equation $ g = f * 1$ is equivalent to $f = g * mu $.

Step by Step Example

This is an example that is computed step by step for those not familiar with the notation used in the definition. Consider the functions $f = mu$ and $g = sigma$. We will now evaluate their convolution $mu * sigma$ at $ n=12$. Their first few terms are listed in the table below.

$$beginarrayccccccccccccc
f & f(1) & f(2) & f(3) & f(4) & f(5) & f(6) & f(7) & f(8) & f(9) & f(10) & f(11) & f(12)
\ hline
mu & 1 & -1 & -1 & 0 & -1 & 1 & -1 & 0 & 0 & 1 & -1 & 0
\
sigma & 1 & 3 & 4 & 7 & 6 & 12 & 8 & 15 & 13 & 18 & 12 & 28
\
endarray$$

The sum iterates over all natural numbers $ d in mathbb N$ that divide $n=12$, thus $d$ assumes all the natural divisors of $n=12 = 2^2cdot 3$. These are $d =1,2,3,4,6,12$. In each summand, we evaluate $g= sigma$ at $d$ and multiply it with $f = mu$ evaluated at $fracnd$. Now we can conclude

$$beginarrayrlccccc
(mu * sigma)(12) &= mu(12)sigma(1) &+mu(6)sigma(2) &+mu(4)sigma(3) &+mu(3)sigma(4) &+mu(2)sigma(6) &+mu(1)sigma(12)
\
&= 0cdot 1 &+ 1cdot 3 &+ 0 cdot 4 &+ (-1)cdot 7 &+ (-1) cdot 12 &+ 1 cdot 28
\
&= 0 & + 3 & 1 0 & -7 & - 12 & + 28
\
&= 12 \
& = id(12)
endarray$$

Lean, 108 100 95 78 75 bytes

def d(f g:_->int)(n):=(list.iota n).foldr(λd s,ite(n%d=0)(s+f d*g(n/d))s)0

Try it online!

More testcases with all of the functions.

Haskell, 46 bytes

(f!g)n=sum[f i*g(div n i)|i<-[1..n],mod n i<1]

Try it online!

Thanks to flawr for -6 bytes and a great challenge! And thanks to H.PWiz for another -6!

Python 3, 59 bytes

lambda f,g,n:sum(f(d)*g(n//d)for d in range(1,n+1)if 1>n%d)

Try it online!

Wolfram Language (Mathematica), 17 bytes

Of course Mathematica has a built-in. It also happens to know many of the example functions. I've included some working examples.

DirichletConvolve

Try it online!