(Description)
令(f(d))表示空间中到原点距离为(d)的整点个数，给定(L,R,k,p)，求
[sum_{d=L}^Rf(d) mathbb{xor} kmod p
]
(L,Rleq 10^{13},R-L+1leq 10^6)。
(Solution)
必然是找规律。OEIS可以直接找到规律，或者打表。
(f)都是(6)的倍数，令(g(x)=frac{f(x)}{6})，可以发现(g(x))为积性函数，即(g(prod p_i^{k^i})=prod g(p_i^{k_i}))。
列出(g(p) （p为质数）)可以知道，若(p\%4=1)，(g(p)=p)；若(p\%4=3)，(g(p)=p+2)。
再看(g(p^2) （p为质数）)可以知道，若(p\%4=1)，(g(p^2)=p^2)；若(p\%4=3)，(g(p^2)=p^2+h(p))，其中(h(p)=8k)，(k)是(p)为第几个(\%4=3)的质数
。
考虑如何更好地表示(h(p))。
由OEIS结论可发现((p-1)h(p)=2(p^2-1))，即(g(p^2)=p^2+frac{2(p^2-1)}{p-1})，且该结论对(p^k)同样成立。
所以
[f(p^k)=egin{cases}1,&p=2\p^k,&p\%4=1\p^k+frac{2(p^k-1)}{p-1},&p\%4=3end{cases}
]
求出(sqrt n)以内的质数，对每个数就可以(O(log n))分解质因数求了。
复杂度(O(Tnlog n))。
PS：注意&的运算优先级比==低。。&运算要加括号！
(ans)求和过程中最多是(10^6)个(2 imes 10^{13})的数，所以还是要取模！
//1261ms 20MB
#include <bits/stdc++.h>
#define pc putchar
#define gc() getchar()
#define pb emplace_back
typedef long long LL;
const int N=3162300,M=1e6+5;
const LL LIM=2e18;
int cnt,P[230000];
bool not_P[N];
inline LL read()
{
LL now=0,f=1; char c=gc();
for(;!isdigit(c);c=='-'&&(f=-1),c=gc());
for(;isdigit(c);now=now*10+c-48,c=gc());
return now*f;
}
void Init(const int n)
{
for(int i=2; i<=n; ++i)
{
if(!not_P[i]) P[++cnt]=i;
for(int j=1; j<=cnt && i*P[j]<=n; ++j)
{
not_P[i*P[j]]=1;
if(!(i%P[j])) break;
}
}
}
void Solve()
{
static LL A[M],f[M];
LL L=read(),R=read(),K=read(),mod=read(),ans=0;
if(!L) ans+=K, ++L;
int n=(int)(R-L+1);
for(int i=1; i<=n; ++i) A[i]=i+L-1,f[i]=1;
for(int i=1,lim=::cnt; i<=lim&&P[i]<=R; ++i)
{
int p=P[i];
LL now=L/p*p; now<L&&(now+=p);
for(now-=L-1; now<=n; now+=p)
{
LL pk=1;
while(!(A[now]%p)) A[now]/=p, pk*=p;
if(p==2) pk=1;
else if((p&3)==3) pk+=2*(pk-1)/(p-1);
f[now]*=pk;
}
}
for(int i=1; i<=n; ++i)
{
if((A[i]&3)==3) f[i]*=(A[i]+2);//本身/剩下为大质数
else f[i]*=A[i];
ans+=((f[i]*6)^K)%mod;
ans>=LIM&&(ans%=mod);//取模！
}
printf("%lld
",ans%mod);
}
int main()
{
Init(N-2);
for(int T=read(); T--; Solve());
return 0;
}