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

新学到了一个求n以内与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}{d}\rfloor} [gcd(\frac{m}{d},i)==1] $

到这一步已经可以递归求了,琪琪说是 \(O(n^{\frac{3}{4}})\) ,不过题解可以继续往下。

④为了方便直接考虑 $ \sum\limits_{i=1}^{n} [gcd(m,i)==1] $ ,反演(大概):

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

⑤交换一下顺序,枚举d,很显然n以内的d的倍数都会贡献一个mu(d):

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


下面的是根据题解的实现,不过是__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;
} else
i3 = 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 Yinku
freopen("Yinku.in", "r", stdin);
#endif // Yinku
sieve();
int T;
cin >> T;
lll n;
while(T--) {
n = read();
cout << S0(n) << endl;
}
}

其实具体的思路还是要先分成两部分来算,但是我当时不会计算这个S1导致T了。其实计算S2的时候在整数较大的时候发生了溢出。也就是c*c的位置。所以说以后非数组的值一律开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;
} else
i3 = 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 Yinku
freopen("Yinku.in", "r", stdin);
#endif // Yinku
sieve();
init();
int T;
cin >> T;
lll n;
while(T--) {
n = read();
cout << S0(n) << endl;
}
}

最新文章

  1. linux中字体的安装以及Terminal字体重叠问题解决
  2. 详解jquery插件中;(function ( $, window, document, undefined )的作用
  3. http协议相关-待续
  4. [.NET领域驱动设计实战系列]专题八:DDD案例:网上书店分布式消息队列和分布式缓存的实现
  5. 字符输出流Writer简要概括
  6. Java多线程与并发库高级应用-Callable与Future的应用
  7. 扩展服务 修改新增Service的默认主题
  8. Android 中 SQLite 数据库的查看
  9. javaScript实现修改输入框之后标红
  10. 黄聪:jquery 校验中国身份证号码
  11. 使用Excel制作万年历(可打印)
  12. linux 常用的基本命令
  13. 【转载】loadrunner使用system()函数调用Tesseract-OCR识别验证码遇到的问题
  14. Android权限安全(1)自定义,检查,使用权限
  15. css透明(支持各浏览器)
  16. WebService只能在本地使用,无法通过网络访问的解决办法
  17. Java基础知识强化之网络编程笔记21:Android网络通信之 Android常用OAuth登录(获取令牌信息)
  18. 让动画不再僵硬:Facebook Rebound Android动画库介绍
  19. vs 中一些快捷键
  20. ORACLE中CHAR、VARCHAR、NVARCHAR

热门文章

  1. 离线下载Express 2015 for Windows 10
  2. MySQL的删除语句
  3. CentOS7 设置电源选项,待机、睡眠、挂起
  4. 1、获取ip地址
  5. scapy - 基于python的数据包操作库
  6. prufer 序列 学习笔记
  7. 域名Whois数据和隐私是最大风险
  8. 【leetcode】1048. Longest String Chain
  9. LeetCode--056--合并区间(java)
  10. 批处理bat文件显示中文乱码解决方式