http://acm.hdu.edu.cn/showproblem.php?pid=6588

新學(xué)到了一個(gè)求n以內(nèi)與m的gcd的和的快速求法。也就是下面的S1。

---

①求：

$ \sum\limits_{i=1}^{n}gcd(m,i) $

②枚舉d：

$ \sum\limits_{d|m} d \sum\limits_{i=1}^{n} [gcd(m,i)==d] $

③顯然：

$ \sum\limits_{d|m} d \sum\limits_{i=1}^{\lfloor\frac{n}ozvdkddzhkzd\rfloor} [gcd(\frac{m}ozvdkddzhkzd,i)==1] $

到這一步已經(jīng)可以遞歸求了，琪琪說是 \(O(n^{\frac{3}{4}})\) ，不過題解可以繼續(xù)往下。

④為了方便直接考慮 $ \sum\limits_{i=1}^{n} [gcd(m,i)==1] $ ，反演（大概）：

$ \sum\limits_{i=1}^{n} \sum\limits_{d|gcd(m,i)} \mu(d) $

⑤交換一下順序，枚舉d，很顯然n以內(nèi)的d的倍數(shù)都會(huì)貢獻(xiàn)一個(gè)mu(d)：

$ \sum\limits_{d|m} \mu(d) \lfloor\frac{n}ozvdkddzhkzd\rfloor $

---

下面的是根據(jù)題解的實(shí)現(xiàn)，不過是__int64的版本。

#include<bits/stdc++.h> using namespace std; typedef long long ll; typedef __int64 lll;const int mod = 998244353; const int MAXN = 10000000;int phi[MAXN + 1]; int pri[MAXN + 1], pritop; bool notpri[MAXN + 1];void sieve() {int n = MAXN;pri[1] = phi[1] = 1;for(int i = 2; i <= n; i++) {if(!pri[i])pri[++pritop] = i, phi[i] = i - 1;for(int j = 1, tmp; j <= pritop && (tmp = i * pri[j]) <= n; j++) {pri[tmp] = 1;if(i % pri[j])phi[tmp] = phi[i] * phi[pri[j]];else {phi[tmp] = phi[i] * pri[j];break;}}} }ll S1(lll n, int m) {//sigma gcd(i,m) [1,n]ll res = 0;for(int T = 1; T * T <= m; ++T) {if(!(m % T)) {res += (n / T) * phi[T];if(T * T != m) {res += (n / (m / T)) * phi[(m / T)];}}}res %= mod;return res; }ll qpow(ll x, int n) {ll res = 1;while(n) {if(n & 1)res = res * x % mod;x = x * x % mod;n >>= 1;}return res; }const int inv2 = qpow(2ll, mod - 2); const int inv6 = qpow(6ll, mod - 2);ll sigma1(ll x) {return x * (x + 1ll) % mod * inv2 % mod; }ll sigma2(ll x) {return x * (x + 1ll) % mod * (2ll * x + 1ll) % mod * inv6 % mod; }ll S2_1(int r, int T) {int c = r / T;ll res = 0;res += 3ll * T * sigma2(c);res += 3ll * sigma1(c);res += c;res %= mod;return res; }ll S2(int r) {ll res = 0;for(int T = 1; T <= r; ++T) {res += 1ll * phi[T] * S2_1(r, T) % mod;}res %= mod;return res; }ll S0(lll n) {lll i, i3;for(i = 1;; ++i) {lll tmp = i * i * i;if(tmp > n) {--i;break;} elsei3 = tmp;}ll res = 0;res += S1(n, i) - S1(i3 - 1, i);res += S2(i - 1);res = (res % mod + mod) % mod;return res; }inline lll read() {lll x = 0;char c;do {c = getchar();} while(c < '0' || c > '9');do {x = (x << 3) + (x << 1) + c - '0';c = getchar();} while(c >= '0' && c <= '9');return x; }int main() { #ifdef Yinkufreopen("Yinku.in", "r", stdin); #endif // Yinkusieve();int T;cin >> T;lll n;while(T--) {n = read();cout << S0(n) << endl;} }

其實(shí)具體的思路還是要先分成兩部分來算，但是我當(dāng)時(shí)不會(huì)計(jì)算這個(gè)S1導(dǎo)致T了。其實(shí)計(jì)算S2的時(shí)候在整數(shù)較大的時(shí)候發(fā)生了溢出。也就是c*c的位置。所以說以后非數(shù)組的值一律開ll就對了。

#include<bits/stdc++.h> using namespace std; typedef long long ll; typedef __int128 lll;const int mod = 998244353; const int MAXN = 10000000;int pk[MAXN + 1]; int sum1[MAXN + 1]; int phi[MAXN + 1]; int pri[MAXN + 1], pritop; bool notpri[MAXN + 1];void sieve() {int n = MAXN;pri[1] = pk[1] = sum1[1] = phi[1] = 1;for(int i = 2; i <= n; i++) {if(!pri[i])pri[++pritop] = i, phi[i] = i - 1, pk[i] = i, sum1[i] = 2 * i - 1;for(int j = 1, p, tmp; j <= pritop && (p = pri[j]) && (tmp = i * p) <= n; j++) {pri[tmp] = 1;if(i % p) {pk[tmp] = pk[p];sum1[tmp] = 1ll * sum1[i] * sum1[p] % mod;phi[tmp] = phi[i] * phi[p];} else {pk[tmp] = pk[i] * p;if(pk[tmp] == tmp) {sum1[tmp] = (1ll * sum1[i] * p % mod + (tmp - tmp / p)) % mod;} else {sum1[tmp] = 1ll * sum1[pk[tmp]] * sum1[tmp / pk[tmp]] % mod;}phi[tmp] = phi[i] * p;break;}}} }int sum2[MAXN + 1];ll qpow(ll x, int n) {ll res = 1;while(n) {if(n & 1)res = res * x % mod;x = x * x % mod;n >>= 1;}return res; }void init() {for(int c = 1, c1 = 2, c2 = 7, f1 = 1; c <= MAXN;) {sum2[c] = ((1ll * c2 - f1 + mod) % mod * sum1[c] % mod * qpow(c, mod - 2) % mod + c + sum2[c - 1]) % mod;++c, ++c1, f1 = c2 + 1;c2 = (1ll * c1 * c1 % mod * c1 % mod - 1 + mod) % mod;} }ll S1(lll n, int m) {//sigma gcd(i,m) [1,n]ll res = 0;for(int T = 1; T * T <= m; ++T) {if(!(m % T)) {res += (n / T) * phi[T];if(T * T != m) {res += (n / (m / T)) * phi[(m / T)];}}}res %= mod;return res; }ll S0(lll n) {lll i, i3;for(i = 1;; ++i) {lll tmp = i * i * i;if(tmp > n) {--i;break;} elsei3 = tmp;}ll res = 0;res += S1(n, i) - S1(i3 - 1, i);res += sum2[i - 1];res = (res % mod + mod) % mod;return res; }inline lll read() {lll x = 0;char c;do {c = getchar();} while(c < '0' || c > '9');do {x = (x << 3) + (x << 1) + c - '0';c = getchar();} while(c >= '0' && c <= '9');return x; }int main() { #ifdef Yinkufreopen("Yinku.in", "r", stdin); #endif // Yinkusieve();init();int T;cin >> T;lll n;while(T--) {n = read();cout << S0(n) << endl;} }

轉(zhuǎn)載于:https://www.cnblogs.com/Yinku/p/11229435.html

總結(jié)

以上是生活随笔為你收集整理的2019 Multi-University Training Contest 1 - 1011 - Function - 数论的全部內(nèi)容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網(wǎng)站內(nèi)容還不錯(cuò)，歡迎將生活随笔推薦給好友。