小A最近开始研究数论题了,这一次他随手写出来一个式子,
∑i=1n∑j=1mgcd(i,j)2,但是他发现他并不太会计算这个式子,你可以告诉他这个结果吗,答案可能会比较大,请模上1000000007。
∑i=1n∑j=1mgcd(i,j)2 ∑d=1min(n,m)d2(∑i=1n∑j=1mgcd(i,j)=d) 记(d)=(∑i=1n∑j=1mgcd(i,j)=d)有F(t)=∑d|tf(d)=[nd][md] 有f(d)=∑d|tμ(td)F(t)=∑d|tμ(td)[nd][md] 处理因子直接计算复杂度 O(nlogn)就可以通过本题 莫比乌斯反演 #include<bits/stdc++.h> using namespace std; typedef long long ll; const int maxn = 1e6; ll prime[maxn + 5]; ll isprime[maxn + 5]; ll mod = 1e9+7; ll mbs[maxn + 5]; int cnt; ll qpow(ll a, ll n, ll mod) { ll ans = 1; a = a % mod; while(n) { if(n&1) { ans = (ans * a) % mod; } a = (a * a) % mod; n >>= 1; } return ans; } void get_mbs(int n) { mbs[1] = 1; isprime[1] = 1; for(int i = 2; i <= n; i++) { if(!isprime[i]) { prime[cnt++] = i; mbs[i] = -1; } for(int j = 0; j < cnt && i * prime[j] <= n; j++) { isprime[i * prime[j]] = 1; if(i % prime[j] == 0) { break; } mbs[i * prime[j]] = -mbs[i]; } mbs[i]+= mbs[i - 1]; } } int main() { ll n, m, k; scanf("%lld%lld", &n, &m); k = 2; get_mbs(maxn); if(n > m) { swap(n, m); } ll ans = 0; for(ll d = 1; d <= n; d++) { ll a = n / d; ll b = m / d; ll res = 0; for(ll i = 1, j = 1; i <= a; i = j + 1) { j = min(a / (a / i), b / (b / i)); ll y = ((a / i) * (b / i)) % mod; res = (res + y * (mbs[j] - mbs[i - 1])) % mod; } ll x = qpow(d, k, mod); ans = (ans + x * res) % mod; } printf("%lld\n", ans); return 0; }