1847 奇怪的數學題

推式子

$$\begin{gathered} \sum _{i=1}^n\sum _{j=1}^{n}sgcd(i,j{)}^k \\ \sum _{d=1}^{n}sgcd(d{)}^k\sum _{i=1}^{\frac{n}{d}}\sum _{j=1}^{\frac{n}{d}}[gcd(i,j)=1] \\ \sum _{d=1}^{n}sgcd(d{)}^k(2\sum _{i=1}^{\frac{n}{d}}?(i)?1) \end{gathered}$$

接下來我們只要求得$sgcd(d{)}^k$的前綴和就可以快樂地數論分塊了，

考慮到這個是求最小質因子，Min_25篩在篩g函數的時候有一步是枚舉最小質因子，

$sgcd=\frac{sgcd}{min(p)}$，對應的就是這一步，所以我們只要舍棄最小質因子，然后統計貢獻就行了，

當然最后我們要加上質數的貢獻，也就是質數的個數。

考慮到要求自然冪數前n項和，這題模數是個頭疼的地方，無法用拉格朗日插值，只能乖乖的寫第二類斯特林數了。

代碼

/*Author : lifehappy */ #pragma GCC optimize(2) #pragma GCC optimize(3) #include <bits/stdc++.h>#define mp make_pair #define pb push_back #define endl '\n' #define mid (l + r >> 1) #define lson rt << 1, l, mid #define rson rt << 1 | 1, mid + 1, r #define ls rt << 1 #define rs rt << 1 | 1using namespace std;typedef long long ll; typedef unsigned long long ull; typedef pair<int, int> pii;const double pi = acos(-1.0); const double eps = 1e-7; const int inf = 0x3f3f3f3f;inline ll read() {ll f = 1, x = 0;char c = getchar();while(c < '0' || c > '9') {if(c == '-') f = -1;c = getchar();}while(c >= '0' && c <= '9') {x = (x << 1) + (x << 3) + (c ^ 48);c = getchar();}return f * x; }#define uint unsigned intconst int N = 1e6 + 10;uint quick_pow(uint a, int n) {uint ans = 1;while(n) {if(n & 1) ans = ans * a;a = a * a;n >>= 1;}return ans; }namespace Sum {uint S[60][60], k;uint Sumk(int n) {uint ret = 0;for(int i = 1; i <= k; ++i) {uint s = 1;for(int j = 0; j <= i; ++j)if((n + 1 - j) % (i + 1)) s *= n + 1 - j;else s *= (n + 1 - j) / (i + 1);ret += S[k][i] *s;}return ret;}void init() {S[0][0] = 1;for(int i = 1; i <= k; i++)for(int j = 1; j <= i; j++)S[i][j] = S[i-1][j-1] + j * S[i-1][j];}}namespace Min_25 {uint prime[N], id1[N], id2[N], m, cnt, k, T;uint a[N], g1[N], sum1[N], g2[N], sum2[N], f[N], ans[N], n;bool st[N];int ID(ll x) {return x <= T ? id1[x] : id2[n / x];}void init() {T = sqrt(n + 0.5);for(int i = 2; i <= T; i++) {if(!st[i]) {prime[++cnt] = i;sum1[cnt] = sum1[cnt - 1] + 1;f[cnt] = quick_pow(i, k);sum2[cnt] = sum2[cnt - 1] + f[cnt];}for(int j = 1; j <= cnt && 1ll * i * prime[j] <= T; j++) {st[i * prime[j]] = 1;if(i % prime[j] == 0) {break;}}}for(ll l = 1, r; l <= n; l = r + 1) {r = n / (n / l);a[++m] = n / l;if(a[m] <= T) id1[a[m]] = m;else id2[n / a[m]] = m;g1[m] = a[m] - 1;g2[m] = Sum::Sumk(a[m]) - 1;}for(int j = 1; j <= cnt; j++) {for(int i = 1; i <= m && 1ll * prime[j] * prime[j] <= a[i]; i++) {g1[i] -= (g1[ID(a[i] / prime[j])] - sum1[j - 1]);g2[i] -= f[j] * (g2[ID(a[i] / prime[j])] - sum2[j - 1]);ans[i] += g2[ID(a[i] / prime[j])] - sum2[j - 1];}}for(int i = 1; i <= m; i++) {ans[i] += g1[i];}}uint solve(int x) {if(x <= 1) return 0;return ans[ID(x)];} }namespace Djs {uint prime[N], phi[N], cnt;bool st[N];void init() {phi[1] = 1;for(int i = 2; i < N; i++) {if(!st[i]) {prime[++cnt] = i;phi[i] = i - 1;}for(int j = 1; j <= cnt && 1ll * i * prime[j] < N; j++) {st[i * prime[j]] = 1;if(i % prime[j] == 0) {phi[i * prime[j]] = phi[i] * prime[j];break;}phi[i * prime[j]] = phi[i] * (prime[j] - 1);}}for(int i = 1; i < N; i++) {phi[i] += phi[i - 1];}}unordered_map<int, uint> ans_s;uint S(int n) {if(n < N) return phi[n];if(ans_s.count(n)) return ans_s[n];uint ans = 1ll * n * (n + 1) / 2;for(uint l = 2, r; l <= n; l = r + 1) {r = n / (n / l);ans -=(r - l + 1) * S(n / l);}return ans_s[n] = ans;} }int main() {// freopen("in.txt", "r", stdin);// freopen("out.txt", "w", stdout);// ios::sync_with_stdio(false), cin.tie(0), cout.tie(0);int n = read(), k = read();Sum::k = k;Sum::init();Djs::init();Min_25::n = n, Min_25::k = k;Min_25::init();uint ans = 0;for(uint l = 1, r; l <= n; l = r + 1) {r = n / (n / l);ans += (Min_25::solve(r) - Min_25::solve(l - 1)) * (2 * Djs::S(n / l) - 1);}cout << ans << endl;return 0; }

總結

以上是生活随笔為你收集整理的1847 奇怪的数学题（杜教筛 + Min_25 + 第二类斯特林数）的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。