狄利克雷卷积

发布时间：2023-05-22

一、介绍

狄利克雷卷积是一种数论函数的一种特殊卷积形式，命名来自于德国数学家彼得·戴利克雷。对于两个数论函数f(n)和g(n)，它们的狄利克雷卷积定义为：

(f * g)(n) = Σd|n f(d)g(n/d)

其中，d是n的正因子。狄利克雷卷积有着广泛的应用，例如在数论中，可以用来定义欧拉函数、除数函数等常见的数论函数；在组合数学中，可以用来计算多项式间的卷积，从而计算组合数的乘积。

二、性质

狄利克雷卷积具有以下的一些性质：

  1. 狄利克雷卷积是交换的，即f * g = g * f。
  2. 对于任意正整数n，(f * g)(n)=0当且仅当f(n')=0或g(n')=0，其中n'为n的任一因子。
  3. 单位函数1与任意函数f卷积的结果仍为f本身，即1 * f = f。
  4. 如果f(1)≠0，那么f与其本身的卷积f * f在n>1的时候恒不为零。

三、应用

1. 求解欧拉函数

欧拉函数是一个关于正整数n的数论函数，表示小于或等于n的正整数中与n互质的数的数目，记为φ(n)。欧拉函数可以使用狄利克雷卷积来定义：

φ(n) = (1 * f)(n)

其中，f(n)表示一个分段函数，在n是合数的情况下，f(n)等于0，在n是质数的情况下，f(n)=n-1。

2. 求解除数函数

除数函数表示一个数的除数个数，记为d(n)。除数函数也可以使用狄利克雷卷积来定义，具体来说，我们需要两个函数：

d(n)= (1 * f)(n)
g(n)=1

其中，f(n)表示一个分段函数，在n是合数的情况下，f(n)等于0，在n是质数的情况下，f(n)=2。 这样，d(n) = (f * g)(n)。

3. 求解多项式卷积

多项式卷积是指给定两个多项式f(x)和g(x)，求解它们的乘积h(x)=f(x)g(x)的过程。这个过程可以使用狄利克雷卷积来实现。 具体来说，我们将f(x)和g(x)中的系数视作数论函数：对于f(x)来说，它的系数函数为f(k)，表示f中次数为k的项的系数；对于g(x)来说，它的系数函数为g(k)，表示g中次数为k的项的系数。这样，h(x)的系数函数就为(f * g)(k)。 有了h(x)的系数函数，我们就可以把它转化为多项式的形式，得到多项式卷积的结果。

四、代码示例

1. 求解欧拉函数

#include <bits/stdc++.h>
using namespace std;
const int N = 1e6;
bool isprime[N+5];
int phi[N+5];
vector<int> primes;
void sieve() {
  memset(isprime, true, sizeof(isprime));
  isprime[1] = false;
  for (int i = 2; i <= N; i++) {
    if (isprime[i]) {
      primes.push_back(i);
      phi[i] = i-1;
    }
    for (auto p : primes) {
      if (p*i > N) break;
      isprime[p*i] = false;
      if (i % p == 0) {
        phi[p*i] = phi[i] * p;
        break;
      } else {
        phi[p*i] = phi[i] * (p-1);
      }
    }
  }
}
int main() {
  sieve();
  for (int i = 1; i <= 20; i++) {
    cout << "phi(" << i << ") = " << phi[i] << endl;
  }
  return 0;
}

2. 求解除数函数

#include <bits/stdc++.h>
using namespace std;
const int N = 1e6;
bool isprime[N+5];
int d[N+5];
vector<int> primes;
void sieve() {
  memset(isprime, true, sizeof(isprime));
  isprime[1] = false;
  d[1] = 1;
  for (int i = 2; i <= N; i++) {
    if (isprime[i]) {
      primes.push_back(i);
      d[i] = 2;
    } else {
      for (auto p : primes) {
        if (p * i > N) break;
        if (i % p == 0) {
          int cnt = 0;
          int x = i;
          while (x % p == 0) {
            cnt++;
            x /= p;
          }
          d[p*i] = d[x] * (cnt+1);
          break;
        } else {
          d[p*i] = d[i] * 2; 
        }
      }
    }
  }
}
int main() {
  sieve();
  for (int i = 1; i <= 20; i++) {
    cout << "d(" << i << ") = " << d[i] << endl;
  }
  return 0;
}

3. 求解多项式卷积

#include <bits/stdc++.h>
using namespace std;
using cd = complex<double>;
const double PI = acos(-1);
const int N = 2e5 + 5;
void fft(vector<cd> &a, bool invert) {
  int n = a.size();
  for (int i = 1, j = 0; i < n; i++) {
    int bit = n >> 1;
    while (j >= bit) {
      j -= bit;
      bit >>= 1;
    }
    j += bit;
    if (i < j) swap(a[i], a[j]);
  }
  for (int len = 2; len <= n; len <<= 1) {
    double ang = 2 * PI / len * (invert ? -1 : 1);
    cd wlen(cos(ang), sin(ang));
    for (int i = 0; i < n; i += len) {
      cd w(1);
      for (int j = 0; j < len/2; j++) {
        cd u = a[i+j], v = w * a[i+j+len/2];
        a[i+j] = u + v;
        a[i+j+len/2] = u - v;
        w *= wlen;
      }
    }
  }
  if (invert) {
    for (int i = 0; i < n; i++) {
      a[i] /= n;
    }
  }
}
vector<int> multiply(vector<int> &a, vector<int> &b) {
  vector<cd> fa(a.begin(), a.end()), fb(b.begin(), b.end());
  int n = 1;
  while (n < a.size() + b.size()) {
    n <<= 1;
  }
  fa.resize(n);
  fb.resize(n);
  fft(fa, false);
  fft(fb, false);
  for (int i = 0; i < n; i++) {
    fa[i] *= fb[i];
  }
  fft(fa, true);
  vector<int> ret(n);
  for (int i = 0; i < n; i++) {
    ret[i] = round(fa[i].real());
  }
  while (!ret.empty() && ret.back() == 0) {
    ret.pop_back();
  }
  return ret;
}
void multiply_inplace(vector<int> &a, vector<int> &b) {
  vector<cd> fa(a.begin(), a.end()), fb(b.begin(), b.end());
  int n = 1;
  while (n < a.size() + b.size()) {
    n <<= 1;
  }
  fa.resize(n);
  fb.resize(n);
  fft(fa, false);
  fft(fb, false);
  for (int i = 0; i < n; i++) {
    fa[i] *= fb[i];
  }
  fft(fa, true);
  a.resize(n);
  for (int i = 0; i < n; i++) {
    a[i] = round(fa[i].real());
  }
  while (!a.empty() && a.back() == 0) {
    a.pop_back();
  }
}
int main() {
  vector<int> a = {1, 2, 3}, b = {4, 5, 6};
  auto c = multiply(a, b);
  for (auto x : c) {
    cout << x << " ";
  }
  cout << endl;
  multiply_inplace(a, b);
  for (auto x : a) {
    cout << x << " ";
  }
  cout << endl;
  return 0;
}
推荐文章