题解：P12828 「DLESS-2」XOR and Number Theory
qczrz6v4nhp6u
·
·
题解
### Solution
考虑一对 $(x,y)$（$x<y$），一个经典结论是 $y-x=(x\oplus y)-2(x\cap\neg y)\le x\oplus y$；而 $\gcd(x,y)=\gcd(y-x,x)\le y-x$，所以 $x\oplus y=\gcd(x,y)$ 当且仅当 $x\oplus y=y-x$ 且 $(y-x)|x$。
不妨设 $d=y-x$ 与 $x=kd,y=(k+1)d$，注意到上面的条件说明了 $x\cap\neg y=0$，所以 $x\subseteq y$，于是有 $d=x\oplus y\subseteq y$。所以说 $(kd,(k+1)d)$ 合法的充要即为 $d\subseteq (k+1)d$。
然后再来看看这个 $x^2\bmod(y^2-xy)$ 是什么东西。直接代入 $x=kd,y=(k+1)d$ 后有：
$$\begin{aligned}
&x^2\bmod(y^2-xy)\\
=&~k^2d^2\bmod((k+1)d^2)\\
=&~d^2(k^2\bmod (k+1))\\
=&~d^2
\end{aligned}$$
那很好了，现在我们要算的答案就是
$$\sum_{d=1}^md^2\sum_{k=2}^{\lfloor\frac nd\rfloor}[d\subseteq kd]$$
注意到 $m$ 很小，考虑直接枚举 $d$，算后面这一部分。
设 $l=\lfloor\log_2d\rfloor+1$，一个显然的观察是 $[d\subseteq kd]$ 是具有至多 $2^l$ 的循环节的，于是我们只需要快速计算 $\sum_{k=1}^r[d\subseteq kd]$（其中 $0\le r< 2^l$）即可得知原式的值。下文为方便先将 $d$ 不断除以 $2$ 直到 $d$ 为奇数，显然不影响答案。
先考虑 $r=2^l-1$（即没有 $r$ 的限制）的情况。考虑从低到高依次确定 $k$ 的第 $i$ 位。若 $d$ 的第 $i$ 位是 $0$，则 $k$ 的这一位显然是 $0$ 是 $1$ 都可以；若 $d$ 的第 $i$ 位是 $1$，则 $k$ 的这一位取决于当前 $kd$ 的第 $i$ 位：若当前 $kd$ 的第 $i$ 位是 $1$，则 $k$ 的第 $i$ 位只能为 $0$，否则只能为 $1$。
直接这样爆搜的复杂度显然为 $O(m^2)$，不过考虑到决策到第 $i$ 位时，只有 $kd$ 的 $i$ 位及以上部分有用，对这个记忆化一下就是 $O(m\sqrt m)$ 了，为啥呢？考虑前 $\frac l2$ 位，只有 $2^{\frac l2}=O(\sqrt d)$ 个状态可以到达；而后 $\frac l2$ 位因为记忆化的缘故也只有 $2^{\frac l2}=O(\sqrt d)$ 个状态，搜索一次的复杂度即为 $O(\sqrt d)$，总复杂度即为 $O(m\sqrt m)$。
上面的讨论都是基于 $r=2^l-1$ 给出的，但是也容易拓展到一般的情况：多设一维 $0/1$ 表示当前 $k\bmod 2^i$ 是否 $\le r\bmod 2^i$ 即可，复杂度仍为 $O(m\sqrt m)$。
### Code
哎哎，怎么跑的又没有官解快。大常数选手落泪了。
```cpp
bool Mst;
#include<bits/stdc++.h>
using namespace std;
using ui=unsigned int;
using ll=long long;
using ull=unsigned long long;
using i128=__int128;
using u128=__uint128_t;
using pii=pair<int,int>;
#define fi first
#define se second
constexpr int N=131075,mod=1e9+7;
inline ll add(ll x,ll y){return (x+=y)>=mod&&(x-=mod),x;}
inline ll Add(ll &x,ll y){return x=add(x,y);}
inline ll sub(ll x,ll y){return (x-=y)<0&&(x+=mod),x;}
inline ll Sub(ll &x,ll y){return x=sub(x,y);}
inline ll qpow(ll a,ll b){
ll res=1;
for(;b;b>>=1,a=a*a%mod)
if(b&1)res=res*a%mod;
return res;
}
int n,m,f[2][N<<1];
inline int dfs(int i,int l,int d,int r,int k,int cur){
int o=(cur>>i)|(1<<(l-i));bool op=(k&((1<<i)-1))<=(r&((1<<i)-1));
if(i==l)return op;
if(~f[op][o])return f[op][o];
int res=0,tmp=(cur+((ll)d<<i))&((1<<l)-1);
if((d>>i&1)<=(cur>>i&1))res+=dfs(i+1,l,d,r,k,cur);
if((d>>i&1)<=(tmp>>i&1))res+=dfs(i+1,l,d,r,k|1<<i,tmp);
return f[op][o]=res;
}
inline void clear(int i,int l,int d,int r,int k,int cur){
int o=(cur>>i)|(1<<(l-i));bool op=(k&((1<<i)-1))<=(r&((1<<i)-1));
if(i==l)return;
if(!~f[op][o])return;
int tmp=(cur+((ll)d<<i))&((1<<l)-1);
if((d>>i&1)<=(cur>>i&1))clear(i+1,l,d,r,k,cur);
if((d>>i&1)<=(tmp>>i&1))clear(i+1,l,d,r,k|1<<i,tmp);
f[op][o]=-1;
}
inline int calc(int d,int r){
if(!r)return 0;
int l=__lg(d)+1;
if(r<(1<<l)){
int res=dfs(0,l,d,r,0,0);
return clear(0,l,d,r,0,0),res;
}
else
return calc(d,(1<<l)-1)*(r/(1<<l))+calc(d,r&((1<<l)-1));
}
void Solve(){
cin>>n>>m;
ll ans=0;
for(int d=1;d<=m;d++){
int x=d;
while(!(x&1))x>>=1;
Add(ans,(ll)d*d%mod*(calc(x,n/d)-1)%mod);
}
cout<<ans<<'\n';
}
bool Med;
int main(){
cerr<<abs(&Mst-&Med)/1048576.0<<endl;
ios::sync_with_stdio(false);
cin.tie(0),cout.tie(0);
memset(f,-1,sizeof f);
int _Test;cin>>_Test;
while(_Test--)Solve();
return 0;
}
```