題目鏈接 "簡單的數學題" 題目描述 輸入一個整數n和一個整數p,你需要求出 $$\sum_{i=1}^n\sum_{j=1}^n (i\cdot j\cdot gcd(i,j))\ mod\ p$$ 其中$gcd(a,b)$表示$a$與$b$的最大公約數 輸入 一行兩個整數$p,n$ 輸出 一行一 ...
題目鏈接
題目描述
輸入一個整數n和一個整數p,你需要求出
\[\sum_{i=1}^n\sum_{j=1}^n (i\cdot j\cdot gcd(i,j))\ mod\ p\]
其中\(gcd(a,b)\)表示\(a\)與\(b\)的最大公約數
輸入
一行兩個整數\(p,n\)
輸出
一行一個整數,為題目中所求值
樣例
樣例輸入
998244353 2000
樣例輸出
883968974
數據範圍
\(n\leq 10^{10}\)
\(5\times 10^8 \leq p \leq 1.1\times 10^9\)
\(p\)為質數(但貌似也可以不是?又不用求逆元)
題解
自己想出來的題!但是連\(WA\)兩發就是因為杜教篩寫掛了……
先不考慮取餘,我們化一下題目中的式子,枚舉\(gcd\)(警告!多公式)。
\[\sum_{i=1}^n\sum_{j=1}^n i\cdot j\cdot gcd(i,j)\]
\[\sum_{d=1}^{n}d\sum_{i=1}^{\left\lfloor \frac{n}{d}\right\rfloor}\sum_{j=1}^{\left\lfloor \frac{n}{d}\right\rfloor}[i\perp j]i\cdot j \cdot d^2\]
\[\sum_{d=1}^{n}d^3\sum_{i=1}^{\left \lfloor \frac{n}{d}\right\rfloor}\sum_{j=1}^{\left\lfloor \frac{n}{d}\right\rfloor}[i\perp j]i\cdot j\]
\[\sum_{d=1}^{n}d^3\sum_{p=1}^{\left \lfloor \frac{n}{d}\right\rfloor}\mu(p)p^2\cdot \Big(\frac{(1+\left\lfloor \frac{n}{dp}\right\rfloor)\left\lfloor \frac{n}{dp}\right\rfloor}{2}\Big)^2\]
額,現在可以使用分塊優化做到\(O(n)\)了,但是這完全不能勝任數據範圍,我們換個角度,設\(dp=T\),枚舉T會有什麼結果?
\[\sum_{T=1}^{n}\Big(\frac{(1+\left\lfloor \frac{n}{T}\right\rfloor)\left\lfloor \frac{n}{T}\right\rfloor}{2}\Big)^2\sum_{d|T}d^3\cdot \mu(\frac{T}{d})(\frac{T}{d})^2\]
\[\sum_{T=1}^{n}\Big(\frac{(1+\left\lfloor \frac{n}{T}\right\rfloor)\left\lfloor \frac{n}{T}\right\rfloor}{2}\Big)^2 T^2\sum_{d|T}d\cdot \mu(\frac{T}{d})\]
現在好像反而變成\(O(n\log n)\)或\(O(n\sqrt{n})\)了,別急,我們看看第二層的求和的意義——狄利克雷捲積,這是\(Id\)函數與\(\mu\)函數的狄利克雷捲積,其值就等於\(\varphi\)。
\[\sum_{T=1}^{n}\Big(\frac{(1+\left\lfloor \frac{n}{T}\right\rfloor)\left\lfloor \frac{n}{T}\right\rfloor}{2}\Big)^2 T^2\varphi(T)\]
現在,我們只需要快速求出一個東西即可——\(T^2\varphi(T)\),前面的部分可以分塊優化,我們急需解決的就是這個函數\(f(T)=T^2\varphi(T)\)的首碼和\(F(T)\)。顯然,這是一個積性函數。
杜教篩的公式:
\[\sum_{i=1}^{n}(f*g)(i)=\sum_{i=1}^{n}\sum_{d|i}f(d)\cdot g(\frac{i}{d})=\sum_{i=1}^{n}g(i)\sum_{j=1}^{\left\lfloor \frac{n}{i}\right\rfloor}f(j)\]
於是我們需要一個函數與\(f\)捲起來,我們根據套路或枚舉發現\(T^2\)項很惱人,於是嘗試把這一項消掉,於是想到了\(g(x)=x^2\)。
\[\sum_{i=1}^{n}\sum_{d|i}d^2\varphi(d)\cdot (\frac{i}{d})^2=\sum_{i=1}^{n}i^2\sum_{j=1}^{\left\lfloor \frac{n}{i}\right\rfloor}f(j)\]
\[\sum_{i=1}^{n}i^2\sum_{d|i}\varphi(d)=\sum_{i=1}^{n}i^2F(\left\lfloor \frac{n}{i}\right\rfloor)\]
根據公式\(\sum_{d|i}\varphi(d)=i\),繼續變形
\[\sum_{i=1}^{n}i^3=F(n)+\sum_{i=2}^{n}i^2F(\left\lfloor \frac{n}{i}\right\rfloor)\]
\[F(n)=\sum_{i=1}^{n}i^3-\sum_{i=2}^{n}i^2F(\left\lfloor \frac{n}{i}\right\rfloor)\]
由於\(p(i)=i^3\)和\(q(i)=i^2\)的首碼和都有公式,我們可以對右邊進行分塊優化,就可以杜教篩了!這道題圓滿解決,時間複雜度\(O(n^{\frac{2}{3}})\)。
不過有些小細節要註意,比如模數乘\(2\)可能會爆\(int\),\(n^2\)可能會爆\(long\ long\),需要先取模再平方
\(Code:\)
#include <map>
#include <cstdio>
#include <cstring>
#include <iostream>
#include <algorithm>
using namespace std;
#define N 5000005
#define ll long long
map<ll, ll>Phi;
ll n, mod, g[N];
int p[N], h[N], phi[N], cnt;
ll sqr(ll x)
{
ll a = 2 * x + 1, b = x + 1, c = x;
if (b % 2 == 0)b /= 2;
else c /= 2;
if (a % 3 == 0)a /= 3;
else
if (b % 3 == 0)b /= 3;
else c /= 3;
a %= mod, b %= mod, c %= mod;
return a * b % mod * c % mod;
}
ll seq(ll x)
{
ll a = x + 1, b = x;
if (a % 2 == 0)a /= 2;
else b /= 2;
a %= mod, b %= mod;
return a * b % mod;
}
ll vas(ll x)
{
ll a = seq(x);
return a * a % mod;
}
ll G(ll x)
{
if (x <= N - 5)
return g[int(x)];
if (Phi.find(x) != Phi.end())
return Phi[x];
ll ans = vas(x);
ll lst = 1;
for (ll i = 2; i <= x; i++)
{
i = x / (x / i);
ll w = (sqr(i) - sqr(lst)) % mod;
ans = (ans - w * G(x / i) % mod) % mod;
lst = i;
}
if (ans < 0)
ans += mod;
Phi.insert(make_pair(x, ans));
return ans;
}
ll Ans(ll x)
{
ll ans = 0, lst = 0;
for (ll i = 1; i <= x; i++)
{
i = x / (x / i);
ll z = seq(x / i);
z = z * z % mod;
ans = (ans + z * (G(i) - G(lst)) % mod) % mod;
lst = i;
}
if (ans < 0)
ans += mod;
return ans;
}
int main()
{
phi[1] = 1;
for (int i = 2; i <= N - 5; i++)
{
if (!h[i])
{
phi[i] = i - 1;
p[++cnt] = i;
}
for (int j = 1; j <= cnt; j++)
{
if (i * p[j] > N - 5)
break;
h[i * p[j]] = 1;
if (i % p[j] == 0)
phi[i * p[j]] = phi[i] * p[j];
else
phi[i * p[j]] = phi[i] * (p[j] - 1);
}
}
cin >> mod >> n;
for (int i = 1; i <= N - 5; i++)
g[i] = (g[i - 1] + 1ll * phi[i] * i % mod * i % mod) % mod;
cout << Ans(n) << '\n';
}