类欧几里得算法小结
基础版
类欧几里得算法被用来处理下面这种式子的值
求值的话我们考虑分类讨论
(1) 若\(a\geq c\)或\(b \geq c\)
注意到\(\lfloor \frac{ac}{b}\rfloor=\lfloor \frac{(a\ mod \ b)*c}{b}\rfloor+c\lfloor \frac{a}{b}\rfloor\),在\(a\geq c\)或\(b\geq c\)时,由于\(\lfloor\frac{a}{c}\rfloor\)或\(\lfloor\frac{b}{c}\rfloor\geq 1\),所以利用上式化简可以得到
(2)若\(a<c\)且\(b<c\)
此时使用上面的式子化简就变得没有意义了,我们需要考虑点其它的东西
考虑我们所计算的式子的实际意义:一条方程为\(ax+cy+b=0\)的直线,在这条直线下方且满足\(x\geq 0\),\(y>0\)的整点的个数,于是我们可以枚举\(y\)来进行化简
记\(m=\lfloor \frac{an+b}{c}\rfloor\)
注意到下取整遇到大于号时可以直接丢掉,化简一下逻辑判断式
带回原式
这样我们就得到了两个十分重要的递推式,注意到它的形式类似于欧几里得算法求\(gcd\),所以被称作类欧几里得算法,时间复杂度和欧几里得算法几乎一致
接下来就是边界条件,这里的边界条件我写的比较多
(1)若\(a=0\),则\(f(a,b,c,n)=(n+1)\lfloor\frac{b}{c}\rfloor\)
(2)若\(n=0\),则\(f(a,b,c,n)=\lfloor\frac{b}{c}\rfloor\)
(3)若\(n=1\),则\(f(a,b,c,n)=\lfloor\frac{a+b}{c}\rfloor+\lfloor\frac{b}{c}\rfloor\)
(4)若\(c<0\),则\(f(a,b,c,n)=f(-a,-b,-c,n)\)
(5)若\(a<0\)或\(b<0\),利用最上面的第一条式子将\(a,b\)取正即可,即\(f(a,b,c,n)=f(a\ mod\ c+c,b\ mod\ c+c,c,n)+\frac{n(n+1)}{2}(\lfloor\frac{a}{c}\rfloor-1)+(n+1)(\lfloor\frac{b}{c}\rfloor-1)\)
反正最后一般就是用到\(a=0\)的情况了
ll work(ll a,ll b,ll c,ll n)
{
if (!a) return (n+1)*(b/c);
if (!n) return b/c;
if (n==1) return (a+b)/c+b/c;
if (c<0) return work(-a,-b,-c,n);
ll d=abs(gcd(gcd(a,b),c));
a/=d;b/=d;c/=d;
if ((a>=c) || (b>=c)) return work(a%c,b%c,c,n)+n*(n+1)/2*(a/c)+(n+1)*(b/c);
if ((a<0) || (b<0)) return work(a%c+c,b%c+c,c,n)+n*(n+1)/2*(a/c-1)+(n+1)*(b/c-1);
ll m=(a*n+b)/c;
return n*m-work(c,c-b-1,a,m-1);
}
扩展版
扩展版主要是两个式子
这两个式子的求解是互相关联的,这一点将在下面的推导过程中体现出来,并且它们的推导和\(f\)有很大的相似之处
求\(g\)
\(g(a,b,c,n)\)较为简单,我们先推导它
(1)若\(a\geq c\)或\(b\geq c\)
和上面一样的化简,注意到此时\(\lfloor\frac{ai}{c}\rfloor\)的前面还有一个系数\(i\),因此按照上面的展开的话得到的系数是一个平方和,同样的\(\lfloor\frac{b}{c}\rfloor\)的系数也由常数\(1\)变成了等差数列求和,这里就直接给出结果
(2) 若\(a<c\)且\(b<c\)
继续将\(\lfloor\frac{ai+b}{c}\rfloor\)展开
为了满足最后的逻辑表达式,\(i\)的取值在\([\lfloor\frac{jc+c-b-1}{a}\rfloor+1,n]\)之中,所以这也是个等差数列,利用求和公式继续化简
好了我们做完了
个鬼啊,那个\(h\)是什么啊
求\(h\)
依然是按照上面的思路分类讨论
(1)若\(a\geq c\)或\(b\geq c\)
由\((a+b+c)^2=a^2+b^2+c^2+2ab+2ac+2bc\)可得
(2)若\(a<c\)且\(b<c\)
首先给出一个变形
于是
好了这回真的做完了
注意到我们需要在\(f,g,h\)之间交替求值,我们对每一对\((a,b,c,n)\)存下它的\(f,g,h\)值即可
node solve(ll a,ll b,ll c,ll n)
{
node ans;
if (!a)
{
ans.f=(n+1)*(b/c)%maxd;
ans.g=(n+1)*n%maxd*inv2%maxd*(b/c)%maxd;
ans.h=(n+1)*(b/c)%maxd*(b/c)%maxd;
return ans;
}
if ((a>=c) || (b>=c))
{
node tmp=solve(a%c,b%c,c,n);
ans.f=(tmp.f+n*(n+1)%maxd*inv2%maxd*(a/c)%maxd+(n+1)*(b/c)%maxd)%maxd;
ans.g=(tmp.g+n*(n+1)%maxd*(n*2+1)%maxd*inv6%maxd*(a/c)%maxd+n*(n+1)%maxd*inv2%maxd*(b/c)%maxd)%maxd;
ans.h=(tmp.h+n*(n+1)%maxd*(n*2+1)%maxd*inv6%maxd*(a/c)%maxd*(a/c)%maxd+(n+1)%maxd*(b/c)%maxd*(b/c)%maxd+tmp.f*(b/c)%maxd*2%maxd+tmp.g*(a/c)%maxd*2%maxd+(a/c)*(b/c)%maxd*n%maxd*(n+1)%maxd)%maxd;
return ans;
}
ll m=(a*n+b)/c;
node tmp=solve(c,c-b-1,a,m-1);m%=maxd;
ans.f=(n*m%maxd+maxd-tmp.f)%maxd;
ans.g=(m*n%maxd*(n+1)%maxd+maxd-tmp.f+maxd-tmp.h)%maxd*inv2%maxd;
ans.h=(n*m%maxd*(m+1)%maxd-tmp.g*2%maxd-tmp.f*2%maxd-ans.f)%maxd;
return ans;
}
记得在最后ans.f=(ans.f+maxd)%maxd;ans.g=(ans.g+maxd)%maxd;ans.h=(ans.h+maxd)%maxd;
强制转为正数即可