传送门
题意:
求 \[ \sum_{i=1}^{n}i^d[gcd(i,n)=1] \]
思路:
我们对上面的式子进行变换, 有:
\[ \begin{aligned} &\sum_{i=1}^{n}i[gcd(i,n)=1]\=&\sum_{i=1}^{n}i\sum_{x|gcd(i,n)}\mu (x)\=&\sum_{i=1}^n i\sum_{x|i,x|n}\mu(x)\=&\sum_{x|n}\mu(x)x^d\sum_{i=1}^{\frac{n}{x}}i^d \end{aligned} \]
以上都是一些套路, 接下来才步入正题.
因为形如 \(\sum_{i=1}^n i^d\) 这种都是一个以 \(n\) 为自变量的, 最高项次数为 \(d+1\) 的多项式.
到这步后, 我们将后面的形式化, 即将多项式表示出来, 设 \(a_i\) 为相关系数, 那么就有式子等于:
\[ \sum_{x|n}\mu(x)x^d\sum_{i=0}^{d+1}a_i\lfloor\frac{n}{x}\rfloor^i \]
变化一下有:
\[ \sum_{i=0}^{d+1}a_i\sum_{x|n}\mu(x)x^d\lfloor\frac{n}{x}\rfloor^i \]
这里面前面的 \(a_i\) 为多项式的系数, 是未知的.
但其实因为我们知道多项式的形式为 \(\sum_{i=0}^{d+1}a_ix^i\), 我们可以直接把 \(x=1,x=2,\cdots,x=d+1\) 的值求出来, 然后高斯消元求解系数.
这里也可以利用拉格朗日插值来求解, 代码是直接抄 https://www.cnblogs.com/cjyyb/p/10503174.html (orz) 的, 思路应该是利用一下等式来求系数:
\[ \sum_{i=0}^n y_i\prod_{i!=j}\frac{x-x_j}{x_i-x_j} \]
那么现在主要就是后面一部分的计算, 我们令 \(f(x)=\mu(x)x^d,g(x)=x^i\), 那么后面一部分可以写为:\(\sum_{x|n}f(x)g(\frac{n}{x})\), 这是狄利克雷卷积的的形式, 因为 \(f,g\) 都为积性, 那么 \(h=f*g(n)\) 也为积性.
所以 \(h(n)=\sum_{x|n}\mu(x)x^d\lfloor\frac{n}{x}\rfloor^i\) 也为积性函数, 那么我们可以考虑单独素因子的贡献. 显然每个素因子只会出现 \(0\) 次或者 \(1\) 次, 否则贡献为 \(0\), 那么有:
\[ \begin{aligned} h(p^a)&=\sum_{j=0}^{a}\mu(p^j)p^{jd}(\frac{p^a}{p^j})^i\&=p^{ai}-p^{ai}p^{d-i} \end{aligned} \]
那么对于每个 \(i,0\leq i\leq d+1\), 求出相应的 \(h(n)\), 再与系数相乘最终结果就出来了.
感觉解法中将多项式形式化出来的想法很巧妙! 没想到多项式还能这么用 hhh, 直接求解多项式的系数也是之前没想到的. 之后对卷积的观察也很重要.
很好的一个题. 细节参考代码:
- /*
- * Author: heyuhhh
- * Created Time: 2019/11/21 19:44:08
- */
- #include <bits/stdc++.h>
- #define MP make_pair
- #define fi first
- #define se second
- #define sz(x) (int)(x).size()
- #define all(x) (x).begin(), (x).end()
- #define INF 0x3f3f3f3f
- #define Local
- #ifdef Local
- #define dbg(args...) do { cout <<#args << "->"; err(args); } while (0)
- void err() { std::cout <<'\n'; }
- template<typename T, typename...Args>
- void err(T a, Args...args) { std::cout <<a << ' '; err(args...); }
- #else
- #define dbg(...)
- #endif
- void pt() {std::cout << '\n'; }
- template<typename T, typename...Args>
- void pt(T a, Args...args) {std::cout <<a << ' '; pt(args...); }
- using namespace std;
- typedef long long ll;
- typedef pair<int, int> pii;
- //head
- const int N = 1000 + 5, MOD = 1e9 + 7;
- ll qpow(ll a, ll b) {
- ll ans = 1;
- while(b) {
- if(b & 1) ans = ans * a % MOD;
- a = a * a % MOD;
- b>>= 1;
- }
- return ans;
- }
- // 求 d 次多项式系数
- struct Lagrange {
- ll f[N], a[N], b[N];
- int d;
- void init(int _d) {
- d = _d;
- //y_i
- for(int i = 1; i <= d + 1; i++) f[i] = (f[i - 1] + qpow(i, d - 1)) % MOD;
- b[0] = 1;
- }
- void work() {
- for(int i = 0; i <= d; i++) {
- for(int j = i + 1; j; j--) b[j] = (b[j - 1] + MOD - 1ll * b[j] * (i + 1) % MOD) % MOD;
- b[0] = 1ll * b[0] * (MOD - i - 1) % MOD;
- }
- for(int i = 0; i <= d; i++) {
- int s = f[i + 1], inv = qpow(i + 1, MOD - 2);
- for(int j = 0; j <= d; j++) if(i != j) s = 1ll * s * qpow((i - j + MOD) % MOD, MOD - 2) % MOD;
- b[0] = 1ll * b[0] * (MOD - inv) % MOD;
- for(int j = 1; j <= d + 1; j++) b[j] = (MOD - 1ll * (b[j] + MOD - b[j - 1]) * inv % MOD) % MOD;
- for(int j = 0; j <= d + 1; j++) a[j] = (a[j] + 1ll * s * b[j]) % MOD;
- for(int j = d + 1; j; j--) b[j] = (b[j - 1] + MOD - 1ll * b[j] * (i + 1) % MOD) % MOD;
- b[0] = 1ll * b[0] * (MOD - i - 1) % MOD;
- }
- }
- }A;
- int d, w;
- ll prod[N];
- void run(){
- A.init(d + 1);
- A.work();
- for(int i = 0; i <= d + 1; i++) prod[i] = 1;
- while(w--) {
- int p, a; cin>> p>> a;
- for(int i = 0; i <= d + 1; i++) {
- ll res = (qpow(p, 1ll * a * i) - qpow(p, 1ll * a * i + d) * qpow(qpow(p, i), MOD - 2) % MOD + MOD) % MOD;
- prod[i] = prod[i] * res % MOD;
- }
- }
- ll ans = 0;
- for(int i = 0; i <= d + 1; i++) ans = (ans + A.a[i] * prod[i] % MOD) % MOD;
- cout <<ans << '\n';
- }
- int main() {
- iOS::sync_with_stdio(false);
- cin.tie(0); cout.tie(0);
- cout << fixed << setprecision(20);
- while(cin>> d>> w) run();
- return 0;
- }
来源: http://www.bubuko.com/infodetail-3298275.html