Dirichlet Convolution

Question

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_\limits{d|n} f\left(\frac{n}{d}\right)\cdot g(d) = \sum_{i\cdot 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 = \frac{n}{d} \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_{n\in\mathbb N}\frac{f(n)}{n^s}\right)\cdot \left(\sum_{n\in\mathbb N}\frac{g(n)}{n^s}\right) = \sum_{n\in\mathbb N}\frac{(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) = \begin{cases}1 & n=1 \\ 0 & n>1 \end{cases}$$ 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) = \begin{cases} (-1)^k & \text{ if } n \text{ is squarefree and } k \text{ is the number of Primefactors of } n \\ 0 & \text{ otherwise } \end{cases} $$ 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) = n\prod_{p|n} \left( 1 - \frac{1}{p}\right) $$ 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 | n} 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 | n} 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) = \begin{cases} 1 & \text{ if } n \text{ is a square number} \\ 0 & \text{otherwise}\end{cases}$$ 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.

$$\begin{array}{c|ccccccccccccc} 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 \\ \end{array}$$

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^2\cdot 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 \$\frac{n}{d}\$. Now we can conclude

$$\begin{array}{rlccccc} (\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) \\ &= 0\cdot 1 &+ 1\cdot 3 &+ 0 \cdot 4 &+ (-1)\cdot 7 &+ (-1) \cdot 12 &+ 1 \cdot 28 \\ &= 0 & + 3 & 1 0 & -7 & - 12 & + 28 \\ &= 12 \\ & = id(12) \end{array}$$

Answer 1

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.

Answer 2

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!

Answer 3

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!

Answer 4

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!

Answer 5

Add++, 51 bytes

D,g,@~,$z€¦~¦*
D,f,@@@,@b[VdF#B]dbRzGb]$dbL$@*z€g¦+

Try it online!

Takes two pre-defined functions as arguments, plus \$n\$, and outputs \$(f * g)(n)\$

How it works

D,g,		; Define a helper function, $g
	@~,	; $g takes a single argument, an array, and splats that array to the stack
		; $g takes the argument e.g. [[τ(x) φ(x)] [3 4]]
		; STACK : 			[[τ(x) φ(x)] [3 4]]
	$z	; Swap and zip:			[[3 τ(x)] [4 φ(x)]]
	€¦~	; Reduce each by execution:	[[τ(3) φ(4)]]
	¦*	; Take the product and return:	τ(3)⋅φ(4) = 4

D,f,		; Define the main function, $f
	@@@,	; $f takes three arguments: φ(x), τ(x) and n (Let n = 12)
		; STACK:			[φ(x) τ(x) 12]
	@	; Reverse the stack:		[12 τ(x) φ(x)]
	b[V	; Pair and save:		[12]			Saved: [τ(x) φ(x)]
	dF#B]	; List of factors:		[[1 2 3 4 6 12]]
	dbR	; Copy and reverse:		[[1 2 3 4 6 12] [12 6 4 3 2 1]]
	z	; Zip together:			[[[1 12] [2 6] [3 4] [4 3] [6 2] [12 1]]]
	Gb]	; Push Saved:			[[[1 12] [2 6] [3 4] [4 3] [6 2] [12 1]] [[τ(x) φ(x)]]]
	$dbL	; Number of dividors:		[[[τ(x) φ(x)]] [[1 12] [2 6] [3 4] [4 3] [6 2] [12 1]] 6]
	$@*	; Repeat:			[[[1 12] [2 6] [3 4] [4 3] [6 2] [12 1]] [[τ(x) φ(x)] [τ(x) φ(x)] [τ(x) φ(x)] [τ(x) φ(x)] [τ(x) φ(x)] [τ(x) φ(x)]]]
	z	; Zip:				[[[τ(x) φ(x)] [1 12]] [[τ(x) φ(x)] [2 6]] [[τ(x) φ(x)] [3 4]] [[τ(x) φ(x)] [4 3]] [[τ(x) φ(x)] [6 2]] [[τ(x) φ(x)] [12 1]]]
	€g	; Run $g over each subarray:	[[4 4 4 6 4 6]]
	¦+	; Take the sum and return:	28

Answer 6

R, 58 bytes

function(n,f,g){for(i in (1:n)[!n%%1:n])F=F+f(i)*g(n/i)
F}

Try it online!

Takes n, f, and g. Luckily the numbers package has quite a few of the functions implemented already.

If vectorized versions were available, which is possible by wrapping each with Vectorize, then the following 45 byte version is possible:

R, 45 bytes

function(n,f,g,x=1:n,i=x[!n%%x])f(i)%*%g(n/i)

Try it online!

Answer 7

APL (Dyalog Classic), 20 bytes

{(⍺⍺¨∘⌽+.×⍵⍵¨)∪⍵∨⍳⍵}

with ⎕IO←1

Try it online!

Easy to solve, hard to test - generally not my type of challenge. Yet, I enjoyed this one very much!

{ } defines a dyadic operator whose operands ⍺⍺ and ⍵⍵ are the two functions being convolved; ⍵ is the numeric argument

∪⍵∨⍳⍵ are the divisors of ⍵ in ascending order, i.e. unique (∪) of the LCMs (∨) of ⍵ with all natural numbers up to it (⍳)

⍵⍵¨ apply the right operand to each

⍺⍺¨∘⌽ apply the left operand to each in reverse

+.× inner product - multiply corresponding elements and sum

---

The same in ngn/apl looks better because of Unicode identifiers, but takes 2 additional bytes because of 1-indexing.

Answer 8

Jelly, 9 bytes

ÆDṚÇ€ḋÑ€Ʋ

Try it online!

Line at the top is the main line of \$f\$, line at the bottom is the main line of \$g\$. \$n\$ is passed as an argument to this function.

Answer 9

Swift 4,  74 70  54 bytes

{n in(1...n).map{n%$0<1 ?f(n/$0)*g($0):0}.reduce(0,+)}

Try it online!

Answer 10

JavaScript (ES6), 47 bytes

Takes input as (f)(g)(n).

f=>g=>h=(n,d=n)=>d&&!(n%d)*f(n/d)*g(d)+h(n,d-1)

Try it online!

Examples

liouville =
n => (-1) ** (D = (n, k = 2) => k > n ? 0 : (n % k ? D(n, k + 1) : 1 + D(n / k, k)))(n)

mobius =
n => (M = (n, k = 1) => n % ++k ? k > n || M(n, k) : n / k % k && -M(n / k, k))(n)

sq =
n => +!((n ** 0.5) % 1)

identity =
n => 1

// sq = liouville * identity
console.log([...Array(25)].map((_, n) => F(liouville)(identity)(n + 1)))

// liouville = mobius * sq
console.log([...Array(20)].map((_, n) => F(mobius)(sq)(n + 1)))

Answer 11

C (gcc), 108 bytes

#define F float
F c(F(*f)(int),F(*g)(int),int n){F s=0;for(int d=0;d++<n;)if(n%d<1)s+=f(n/d)*g(d);return s;}

Straightforward implementation, shamelessly stolen from Leaky Nun's Python answer.

Ungolfed:

float c(float (*f)(int), float (*g)(int), int n) {
    float s = 0;
    for(int d = 1; d <= n;++d) {
        if(n % d == 0) {
            s += f(n / d) * g(d);
        }
    }
    return s;
}

Try it online!

Answer 12

F#, 72 bytes

let x f g n=Seq.filter(fun d->n%d=0){1..n}|>Seq.sumBy(fun d->f(n/d)*g d)

Takes the two functions f and g and a natural number n. Filters out the values of d that do not naturally divide into n. Then evaluates f(n/d) and g(d), multiples them together, and sums the results.

Answer 13

Pari/GP, 32 bytes

(f,g,n)->sumdiv(n,d,f(n/d)*g(d))

Try it online!

There is a built-in dirmul function, but it only supports finite sequences.

Answer 14

APL(NARS), 47 chars, 94 bytes

{(m⍵[1])×n⍵[2]}{+/⍺⍺¨{k←⍳⍵⋄(⍵÷b),¨b←k/⍨0=k∣⍵}⍵}

where m and n are the function one has to use (this because i don't know how to call one array of function in one function in APL). Using the example above the multiplication of Mobius function (here it is 12π) and sum of divisors function (here it is 11π) for value 12 that multiplication would be:

  {(12π⍵[1])×11π⍵[2]}{+/⍺⍺¨{k←⍳⍵⋄(⍵÷b),¨b←k/⍨0=k∣⍵}⍵}12
12

if one has to calculate some other value:

  {(12π⍵[1])×11π⍵[2]}{+/⍺⍺¨{k←⍳⍵⋄(⍵÷b),¨b←k/⍨0=k∣⍵}⍵}1002
1002
  {(12π⍵[1])×11π⍵[2]}{+/⍺⍺¨{k←⍳⍵⋄(⍵÷b),¨b←k/⍨0=k∣⍵}⍵}1001
1001
  {(12π⍵[1])×11π⍵[2]}{+/⍺⍺¨{k←⍳⍵⋄(⍵÷b),¨b←k/⍨0=k∣⍵}⍵}20000x
20000 

One can see if for example the first 2000 number the function result is the identity

  (⍳2000)≡{(12π⍵[1])×11π⍵[2]}{+/⍺⍺¨{k←⍳⍵⋄(⍵÷b),¨b←k/⍨0=k∣⍵}⍵}¨⍳2000
1