BZOJ3601 一個人的數論 【數論 + 高斯消元】
阿新 • • 發佈：2018-06-25
HA urn 個人 mat CM ali play clu ons
我們記
\[g(x) = \sum\limits_{i = 1}^{x}i^{d}\]
然後就是最匪夷所思的地方，我們大力猜想這是關於\(x\)的一個\(d + 1\)次多項式
即
\[g(x) = \sum\limits_{i = 1}^{d + 1}a_ix^{i}\]
只需高斯消元得出系數\(a_i\)
後面是一個狄利克雷卷積
\(F(x) = \mu(x)x^{d}\)是一個積性函數，\(F(x) = x^{i}\)顯然也是一個積性函數
兩個積性函數的狄利克雷卷積依舊是一個積性函數
所以我們只需計算\(n\)的所有質因子的函數值乘起來
所以我們記
\[h(p^{k}) = \sum\limits_{c|p^{k}}\mu(c)c^{d}(\frac{p^{k}}{c})^{i}\]
顯然只有\(\mu(1)\)和\(\mu(p)\)兩項
化簡得
\[h(p^{k}) = p_{ki}(1 - p^{d - i})\]
可以\(O(1)\)計算
\(O(dw)\)計算即可
題目鏈接
BZOJ3601
題解
挺神的
首先有
\[
\begin{aligned}
f(n) &= \sum\limits_{x = 1}^{n} x^{d} [(x,n) = 1] \&= \sum\limits_{x = 1}^{n} x^{d} \sum\limits_{c|(x,n)}\mu(c) \&= \sum\limits_{c|n}\sum\limits_{x = 1}^{\frac{n}{c}} (cx)^{d} \mu(c) \&= \sum\limits_{c|n}\mu(c)c^{d}\sum\limits_{x = 1}^{\frac{n}{c}} x^{d} \\end{aligned}
\]
我們記
\[g(x) = \sum\limits_{i = 1}^{x}i^{d}\]
然後就是最匪夷所思的地方，我們大力猜想這是關於\(x\)的一個\(d + 1\)次多項式
即
\[g(x) = \sum\limits_{i = 1}^{d + 1}a_ix^{i}\]
只需高斯消元得出系數\(a_i\)
然後\(f(n)\)可以繼續化簡
\[
\begin{aligned}
f(n) &= \sum\limits_{c|n}\mu(c)c^{d}g(\frac{n}{c}) \&= \sum\limits_{c|n}\mu(c)c^{d}\sum\limits_{i = 1}^{d + 1} a_i(\frac{n}{c})^{i} \&= \sum\limits_{i = 1}^{d + 1}a_i\sum\limits_{c|n}\mu(c)c^{d}(\frac{n}{c})^{i}
\end{aligned}
\]
後面是一個狄利克雷卷積
\(F(x) = \mu(x)x^{d}\)是一個積性函數，\(F(x) = x^{i}\)顯然也是一個積性函數
兩個積性函數的狄利克雷卷積依舊是一個積性函數
所以我們只需計算\(n\)的所有質因子的函數值乘起來
所以我們記
\[h(p^{k}) = \sum\limits_{c|p^{k}}\mu(c)c^{d}(\frac{p^{k}}{c})^{i}\]
顯然只有\(\mu(1)\)和\(\mu(p)\)兩項
化簡得
\[h(p^{k}) = p_{ki}(1 - p^{d - i})\]
可以\(O(1)\)計算
所以式子就化為
\[f(n) = \sum\limits_{i = 1}^{d + 1}a_i\prod_{i=1}^{w}h(p_i^{k_i})\]
\(O(dw)\)計算即可
總復雜度\(O(d^3 + dw)\)
#include<algorithm>
#include<iostream>
#include<cstring>
#include<cstdio>
#include<cmath>
#include<map>
#define Redge(u) for (int k = h[u],to; k; k = ed[k].nxt)
#define REP(i,n) for (int i = 1; i <= (n); i++)
#define mp(a,b) make_pair<int,int>(a,b)
#define cls(s) memset(s,0,sizeof(s))
#define cp pair<int,int>
#define LL long long int
using namespace std;
const int maxn = 105,maxm = 1005,INF = 1000000000,P = 1000000007;
inline int read(){
int out = 0,flag = 1; char c = getchar();
while (c < 48 || c > 57){if (c == ‘-‘) flag = -1; c = getchar();}
while (c >= 48 && c <= 57){out = (out << 3) + (out << 1) + c - 48; c = getchar();}
return out * flag;
}
int w,d,p[maxm],k[maxm],a[maxn];
int A[maxn][maxn],N;
inline int qpow(int a,LL b){
if (b < 0) b += P - 1;
int re = 1;
for (; b; b >>= 1,a = 1ll * a * a % P)
if (b & 1) re = 1ll * re * a % P;
return re;
}
void gause(){
for (int i = 1; i <= N; i++){
int j = i;
/*for (int k = i + 1; k <= N; k++)
if (A[k][i] > A[j][i]) j = k;
if (j != i) for (int k = i; k <= N + 1; k++) swap(A[j][k],A[i][k]);*/
for (j = i + 1; j <= N; j++){
int t = 1ll * A[j][i] * qpow(A[i][i],P - 2) % P;
for (int k = i; k <= N + 1; k++)
A[j][k] = ((A[j][k] - 1ll * A[i][k] * t % P) % P + P) % P;
}
}
for (int i = N; i; i--){
for (int j = i + 1; j <= N; j++)
A[i][N + 1] = ((A[i][N + 1] - 1ll * a[j] * A[i][j] % P) % P + P) % P;
a[i] = 1ll * A[i][N + 1] * qpow(A[i][i],P - 2) % P;
}
}
void cal(){
N = d + 1;
for (int x = 1; x <= N; x++){
A[x][N + 1] = (A[x - 1][N + 1] + qpow(x,d)) % P;
for (int j = 1; j <= N; j++) A[x][j] = qpow(x,j);
}
gause();
int s1 = 0,s2 = 0;
for (int i = 1; i <= N; i++) s1 = (s1 + 1ll * a[i] * qpow(5,i) % P) % P;
for (int i = 1; i <= 5; i++) s2 = (s2 + qpow(i,d)) % P;
}
void work(){
int ans = 0;
for (int i = 1; i <= N; i++){
int tmp = a[i];
for (int j = 1; j <= w; j++)
tmp = 1ll * tmp * qpow(p[j],1ll * k[j] * i) % P * (((1 - qpow(p[j],d - i)) % P + P) % P) % P;
ans = (ans + tmp) % P;
}
printf("%d\n",ans);
}
int main(){
d = read(); w = read();
REP(i,w) p[i] = read(),k[i] = read();
cal();
work();
return 0;
}BZOJ3601 一個人的數論 【數論 + 高斯消元】