简化题面
求
[sum_{i=1}^nsum_{j=1}^mvarphi{(i imes j)}
]
答案对(998244353)取模。多组询问
(1le Tle 10^4)，(1le n，mle 10^5)
解题思路
这题在我计划里面躺了最少有五个月了，当时觉得好难啊就没做= =，最近拿出来看看好简单啊于是就做了
考虑大力推柿子
对于(varphi{(i imes j)})，我们有一个不是很显然的也不是很常见的结论，
[varphi{(i imes j)}=frac{varphi{(i)} imes varphi{(j)} imes gcd(i,j)}{varphi{(gcd(i,j))}}
]
考虑如何证明，我们将(varphi)都用定义拆开，可以写成比较鬼畜的形式
[i imes jprod_{p|i or j} frac{p-1}{p}=frac{iprod_{p|i} frac{p-1}{p} jprod_{p|j} frac{p-1}{p}gcd(i,j)}{gcd(i,j)prod_{p|i and j} frac{p-1}{p}}
]
显然两边是相等的。
那么既然已经知道了这玩意那我们直接就可以做了
[sum_{i=1}^nsum_{j=1}^m frac{varphi{(i)} imes varphi{(j)} imes gcd(i,j)}{varphi{(gcd(i,j))}}
]
[=sum_{p=1}^n frac{p}{varphi{(p)}}sum_{i=1}^nsum_{j=1}^m varphi{(i)}varphi{(j)} [gcd(i,j)=p]
]
[=sum_{p=1}^n frac{p}{varphi{(p)}} sum_{i=1}^{frac{n}{p}}sum_{j=1}^{frac{m}{p}} varphi{(ip)}varphi{(jp)} [gcd(i,j)=1]
]
[=sum_{p=1}^n frac{p}{varphi{(p)}} sum_{i=1}^{frac{n}{p}}sum_{j=1}^{frac{m}{p}} varphi{(ip)}varphi{(jp)} sum_{d|i,j}mu{(d)}
]
[=sum_{p=1}^n frac{p}{varphi{(p)}} sum_{d=1}^n mu{(d)} sum_{i=1}^{frac{n}{dp}}sum_{j=1}^{frac{m}{dp}} varphi{(idp)}varphi{(jdp)}
]
[=sum_{Q=1}^n sum_{p|Q} frac{pmu{(frac{Q}{p})}}{varphi{(p)}} sum_{i=1}^{frac{n}{Q}} varphi{(iQ)}sum_{j=1}^{frac{m}{Q}}varphi{(jQ)}
]
显然到这里没办法继续推了，考虑如何进行预处理。
第一部分显然可以在(O(nlogn))的时间内进行预处理得到
那么考虑如何处理后面的玩意，我们令(g(k,n)=sum_{i=1}^n varphi{(i imes k)})，这个玩意显然是可以递推的(g(i,j)=g(i,j-1)+varphi{(i imes j)})
于是这个玩意又可以在(O(nlogn))的时间复杂度内预处理出来
那么原柿子现在变成了
[sum_{Q=1}^n sum_{p|Q} frac{pmu{(frac{Q}{p})}}{varphi{(p)}} g(Q,frac{n}{Q}) g(Q,frac{m}{Q})
]
那么到这里一个很自然的想法是我们可以进行分块打表预处理
考虑设(t(i,j,Q))，前面的(i,j)是分别处于的块(frac{n}{Q})和(frac{m}{Q})，这个玩意还是可以递推的(t(i,j,Q)=t(i,j,Q-1)+sum_{p|Q} frac{pmu{(frac{Q}{p})}}{varphi{(p)}} g(Q,i) g(Q,j))
我们直接设定大小为(len)的块暴力(O(n^2))预处理就可以了，最后统计答案的时候散块暴力整块直接分块统计答案就可以了。
时间复杂度是(O(nlogn+nlen^2+T(sqrt n+frac{n}{len})))
(mathcal{Code})
// Author: Ame__
#include<bits/stdc++.h>
#include<stdint.h>
#define _ 0
#define AME__DEBUG
#define bomb exit(0)
#define LOG(FMT...) fprintf(stderr , FMT)
#define TOWA(FMT...) fprintf(stdout , FMT)
using namespace std;
/*Grievous Lady*/
typedef int32_t i32;
typedef int64_t i64;
typedef double qwq;
const int BUF_SIZE = 1 << 12;
char buf[BUF_SIZE] , *buf_s = buf , *buf_t = buf + 1;
#define PTR_NEXT()
{
buf_s ++;
if(buf_s == buf_t)
{
buf_s = buf;
buf_t = buf + fread(buf , 1 , BUF_SIZE , stdin);
}
}
#define mians(_s_)
{
while(!isgraph(*buf_s)) PTR_NEXT();
char register *_ptr_ = (_s_);
while(isgraph(*buf_s) || *buf_s == '-')
{
*(_ptr_ ++) = *buf_s;
PTR_NEXT();
}
(*_ptr_) = ' ';
}
template <typename _n_> void mian(_n_ & _x_){
while(*buf_s != '-' && !isdigit(*buf_s)) PTR_NEXT();
bool register _nega_ = false; if(*buf_s == '-'){ _nega_ = true; PTR_NEXT(); }
_x_ = 0; while(isdigit(*buf_s)){ _x_ = _x_ * 10 + *buf_s - '0'; PTR_NEXT(); } if(_nega_) _x_ = -_x_;
}
const i32 kato = 1e5 + 10;
const i32 atri = 1e5;
const i32 deco = 1e2 + 10;
const i32 mod = 998244353;
template <typename _n_> bool cmax(_n_ &a , const _n_ &b){ return a < b ? a = b , 1 : 0; }
template <typename _n_> bool cmin(_n_ &a , const _n_ &b){ return a > b ? a = b , 1 : 0; }
i32 prime[kato] , phi[kato] , mu[kato] , inv[kato] , sum[kato];
bool ispri[kato];
i32 n , m , cnt , T , b = 50;
i32 *g[kato] , *t[deco][deco];
inline i32 quick_pow(i32 a , i32 b){
i32 res = 1;
for(; b ; b >>= 1 , a = static_cast<i64>(a) * a % mod){
if(b & 1){
res = static_cast<i64>(res) * a % mod;
}
}
return res;
}
inline void pri(){
for(i32 i = 2;i <= atri;i ++){
if(!ispri[i]){
prime[++ cnt] = i;
phi[i] = i - 1;
mu[i] = -1;
}
for(i32 j = 1;j <= cnt && i * prime[j] <= atri;j ++){
ispri[i * prime[j]] = 1;
if(i % prime[j] == 0){
phi[i * prime[j]] = phi[i] * prime[j];
break;
}
phi[i * prime[j]] = phi[i] * phi[prime[j]];
mu[i * prime[j]] = -mu[i];
}
inv[i] = quick_pow(phi[i] , mod - 2);
}
for(i32 i = 1;i <= atri;i ++){
for(i32 j = i , d = 1;j <= atri;j += i , d ++) sum[j] = ((sum[j] + static_cast<i64>(i) * inv[i] % mod * mu[d] % mod) + mod) % mod;
}
for(i32 i = 1;i <= atri;i ++){
g[i] = new i32[atri / i + 10] , g[i][0] = 0;
for(i32 j = 1;j <= atri / i;j ++) g[i][j] = (g[i][j - 1] + phi[i * j]) % mod;
}
for(i32 i = 1;i <= b;i ++){
for(i32 j = i;j <= b;j ++){
i32 len = atri / max(i , j);
t[i][j] = new i32[len + 10] , t[i][j][0] = 0;
for(i32 k = 1;k <= len;k ++) t[i][j][k] = (t[i][j][k - 1] + static_cast<i64>(sum[k]) * g[k][i] % mod * g[k][j] % mod) % mod;
}
}
}
inline int Ame_(){
#ifdef AME__
freopen(".in" , "r" , stdin); freopen(".out" , "w" , stdout); int nol_cl = clock();
#endif
mian(T); phi[1] = mu[1] = inv[1] = 1; pri();
for(; T --> 0 ;){
mian(n) , mian(m);
if(n > m) swap(n , m);
i32 ans = 0;
for(i32 i = 1;i <= m / b;i ++) ans = (ans + static_cast<i64>(sum[i]) * g[i][n / i] % mod * g[i][m / i] % mod) % mod;
for(i32 l = m / b + 1 , r;l <= n;l = r + 1){
r = min(n / (n / l) , m / (m / l));
ans = ((ans + t[n / l][m / l][r] - t[n / l][m / l][l - 1]) % mod + mod) % mod;
}
TOWA("%d
" , ans);
}
#ifdef AME__TIME
LOG("Time: %dms
", int((clock() - nol_cl) / (qwq)CLOCKS_PER_SEC * 1000));
#endif
return ~~(0^_^0); /*さようならプログラム*/
}
int Ame__ = Ame_();
int main(){;}