【BZOJ3817/UOJ42】Sum（类欧）

题面

题解

令$x=\sqrt r$，那么要求的式子是$$\sum_{d=1}^n(-1){[dx]}$$

不难发现，对于每个$d$而言的取值只和$[dx]$的奇偶性相关。

如果$x$是个整数，也就是$r$是完全平方数的时候，显然是可以直接算答案的。

计算答案的时候显然之和有几个奇数或者几个偶数相关（只要求一个另外一个就是补集）

比如说我们来求有几个是偶数，那么要满足的条件就是$[dx]=2*[\frac{dx}{2}]$。

先重新写下式子，我们写成这个样子$$\sum_{d=1}^n(1-2*([dx]\ mod\ 2)$$

这个显然成立，当$[dx]$是偶数的时候贡献是$1$，奇数的时候贡献是$0$。

而$[dx]\ mod\ 2$可以写成减法的形式。所以原式可以写成

\[\sum_{d=1}^n(1-2*([dx]-2*[\frac{dx}{2}]))
\]

化简之后得到了

\[n+\sum_{d=1}^n (4*[\frac{dx}{2}]-2*[dx])
\]

现在把模型转化一下，变得更加一般一点，那么要求解的东西是$$\sum _{d=1}^n[\frac{ax+b}{c}d]$$

听说这个玩意叫做类欧？类欧几里得算法。

令$k=\frac{ax+b}{c}$，那么式子可以化简为$\sum[kd]$

这里进行分类讨论

1.当$k>=1$的时候

\[\begin{aligned}\sum_{d=1}^n[kd]&=[d*(k-[\frac{ax+b}{c}])]+d*[\frac{ax+b}{c}]\end{aligned}
\]

显然后面那一半是可以直接计算的，而前面这一半令$k'=k-[\frac{ax+b}{c}]$

我们可以递归处理。

2.$k<1$

转化到一个坐标系上面去，我们要求的东西本质上就是有一条直线$y=kx$，要求解$x\in[1,n]$时，与$x=n$、$x$正半轴围成的三角形内部整点的个数。

我们把这个三角形逆时针旋转$90$度，再沿着$y$轴翻转过来，让长度为$n$的那条边靠着$y$轴，这样子翻转过来之后，$k$变成了倒数，$n$变成了$[nk]$，然后拿矩形减去多出来的部分，不就是矩形减去$k>=1$的情况了吗？

再讨论一下取倒数的结果，斜率$k=\frac{ax+b}{c}$，取倒数之后$k=\frac{c}{ax+b}$，然后分母有理化一下$k=\frac{c(ax-b)}{a^2r-b^2}$。

总的时间复杂度是$log$级别的。

#include<iostream>

#include<cstdio>

#include<algorithm>

using namespace std;

#define ll long long

inline int read()

{

	int x=0;bool t=false;char ch=getchar();

	while((ch<'0'||ch>'9')&&ch!='-')ch=getchar();

	if(ch=='-')t=true,ch=getchar();

	while(ch<='9'&&ch>='0')x=x*10+ch-48,ch=getchar();

	return t?-x:x;

}

double x;ll n,r;

ll Calc(ll a,ll b,ll c,ll n)

{

	if(!n)return 0;

	ll d=__gcd(a,__gcd(b,c));a/=d;b/=d;c/=d;

	ll k=(a*x+b)/c;

	if(!k)

	{

		ll N=(a*x+b)/c*n;

		return n*N-Calc(a*c,-b*c,a*a*r-b*b,N);

	}

	else

		return k*n*(n+1)/2+Calc(a,b-c*k,c,n);

}

int main()

{

	int T=read();

	while(T--)

	{

		n=read();r=read();x=sqrt(r);

		ll k=x;

		if(k*k==r)

		{

			if(k&1)printf("%lld\n",n-2*((n+1)/2));

			else printf("%lld\n",n);

		}

		else printf("%lld\n",n+Calc(1,0,2,n)*4-Calc(1,0,1,n)*2);

	}

	return 0;

}