HDU 6706 杜教筛

传送门: http://acm.hdu.edu.cn/showproblem.php?pid=6706

题意

$$求解\sum_{i=1}^n\sum_{j=1}^igcd(i^a-j^a,i^b-j^b)[gcd(i,j)=1]$$

思路

$$gcd(i^a-j^a,i^b-j^b)=i^{gcd(a,b)}-j^{gcd(a,b)}$$

题目说gcd(a,b)=1,所以那个复杂的式子直接变成i-j.

$$\sum_{i=1}^n\sum_{j=1}^i(i-j)[gcd(i,j)=1]$$

$$\sum_{i=1}^n\sum_{j=1}^ii[gcd(i,j)=1]-\sum_{i=1}^n\sum_{j=1}^ij[gcd(i,j)=1]$$

左边为与i互质的个数,右边为比i小并且互质的数字和,即为:

$$\sum_{i=1}^ni\varphi (i)-\frac{\sum_{i=1}^ni\varphi (i)+1}{2}$$

$$\frac{\sum_{i=1}^ni\varphi (i)-1}{2}$$

因为n有1e9,所以需要杜教筛优化。

$$f*g(n)=\sum_{dn}f(d)g(\frac{n}{d})=\sum_{dn}id(d)*\varphi (d)*g(\frac{n}{d})$$

很显然,设g=id,则:

$$f*g(n)=\sum_{dn}\varphi (d)id(n)=n\sum_{dn}\varphi (d)=id(n)(\varphi *I)=id^2(n)=n^2$$

$$所以得:S(n)=\sum_{i=1}^n(f*g)i-\sum_{i=2}^ng(i)S(\frac{n}{i})=\frac{n(n+1)(2n+1)}{6}-\sum_{i=2}^niS(\frac{n}{i})$$

Code

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
##include "bits/stdc++.h"
using namespace std;

typedef long long ll;
typedef pair<int, int> pii;
typedef pair<double, double> pdd;
typedef pair<ll, ll> pll;

##define endl "\n"
##define eb emplace\_back
##define mem(a, b) memset(a , b , sizeof(a))

const ll INF = 0x3f3f3f3f;
// const ll mod = 998244353;
const ll mod = 1e9 + 7;
const double eps = 1e-6;
const double PI = acos(-1);
const double R = 0.57721566490153286060651209;

template<typename T>
inline void read(T &a) {char c = getchar();T x = 0, f = 1;
while (!isdigit(c)) {if (c == '-')f = -1;c = getchar();}
while (isdigit(c)) {x = (x << 1) + (x << 3) + c - '0';c = getchar();}a = f * x;
}

const int N = 1e6 + 10;

bool is\_prime[N];
int prime[N], cnt;
int phi[N];
ll sum\_i\_phi[N];

void init() {
phi[1] = 1;
for(int i = 2;i < N; i++) {
if(!is\_prime[i]) prime[++cnt] = i, phi[i] = i - 1;
for(int j = 1;j <= cnt && i * prime[j] < N; j++) {
is\_prime[i * prime[j]] = 1;
if(i % prime[j] == 0) {
phi[i * prime[j]] = phi[i] * prime[j];
break;
}
else phi[i * prime[j]] = phi[i] * phi[prime[j]];
}
}
for(int i = 1;i < N; i++) {
sum\_i\_phi[i] = (sum\_i\_phi[i - 1] + 1ll * i * phi[i] % mod) % mod;
}
}

ll quick\_pow(ll a, ll b) {
ll ans = 1;
while(b) {
if(b & 1) ans = ans * a % mod;
a = a * a % mod;
b >>= 1;
}
return ans % mod;
}

ll inv2;
ll inv6;

map<ll, ll> mp;

ll S(ll x) {
if(x < N) return sum\_i\_phi[x];
else if(mp[x]) return mp[x];
ll ans = x * (x + 1) % mod * (2 * x + 1) % mod * inv6 % mod;
for(int l = 2, r;l <= x; l = r + 1) {
r = min(x, x / (x / l));
ans = (ans - 1ll * (r + l) % mod * (r - l + 1) % mod * inv2 % mod * S(x / l)) % mod;
}
return mp[x] = ans;
}

void solve() {
init();
inv2 = quick\_pow(2, mod - 2);
inv6 = quick\_pow(6, mod - 2);
int \_; scanf("%d",&\_);
while(\_--) {
ll n, a, b; scanf("%lld%lld%lld",&n,&a,&b);
ll ans = (S(n) - 1) % mod;
ans = ans * inv2 % mod;
printf("%lld\n",(ans % mod + mod) % mod);
}
}

signed main() {
ios\_base::sync\_with\_stdio(false);
// cin.tie(nullptr);
// cout.tie(nullptr);
##ifdef FZT\_ACM\_LOCAL
freopen("in.txt", "r", stdin);
freopen("out.txt", "w", stdout);
signed test\_index\_for\_debug = 1;
char acm\_local\_for\_debug = 0;
do {
if (acm\_local\_for\_debug == '$') exit(0);
if (test\_index\_for\_debug > 20)
throw runtime\_error("Check the stdin!!!");
auto start\_clock\_for\_debug = clock();
solve();
auto end\_clock\_for\_debug = clock();
cout << "Test " << test\_index\_for\_debug << " successful" << endl;
cerr << "Test " << test\_index\_for\_debug++ << " Run Time: "
<< double(end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s" << endl;
cout << "--------------------------------------------------" << endl;
} while (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug));
##else
solve();
##endif
return 0;
}

本文作者:jujimeizuo
本文地址https://blog.jujimeizuo.cn/2021/03/09/hdu-6706/
本博客所有文章除特别声明外,均采用 CC BY-SA 3.0 协议。转载请注明出处!