题解 P5170 【【模板】类欧几里得算法】
ix35
·
·
题解
最近在学万能欧几里得的套路,就记录一下。
目前的十篇题解中只有一篇是用万能欧几里得的,我感觉讲的也不是非常清楚,就自己写一遍吧。
P5170 【模板】类欧几里得算法
解决什么?
很多类欧几里得问题都可以归结为以下这个模型:
考虑笛卡尔坐标系中一条直线 y=\dfrac{Px+R}{Q},并且做出所有形如 x=k\ (k\in Z) 的竖线以及 y=h\ (h\in Z) 的横线。
假想从左往右扫描这条直线,每遇到这条直线与一个横线的交点,则执行 U 操作;每遇到这条直线与一个竖线的交点,则执行 R 操作,规定同时遇到一条横线和竖线(即到达一个整点)时先执行 U 后执行 R。
这样一来,我们可以将一条直线转译成一个由 U,R 组成的操作序列。
例如,考虑 y=\dfrac{3x+2}{5},其图像如下:
在 x\in (0,6] 的范围内,其对应的操作序列就是:URRURRURUR。
同时,操作要满足可合并性。例如上面的操作序列,我们设 P=UR,可以进一步转化为 PRPRPP。
操作的内容可能是对于一些变量做修改,具体的情况我们最后再说。而万能欧几里得的套路可以帮我们快速计算进行这样的操作序列后的某些信息。
我们将两个操作序列 S 和 T 的合并序列称为 S+T 或 ST,这个合并过程不一定满足交换律。
如何解决?
下面这段里面有作为整数的 R 和作为操作的 R,但好像并不会出现什么混乱,这里就不改了。
方便起见,我们用一个“结点”(Node)来刻画一个操作序列。
解决问题的过程是递归的,我们设 solve(P,Q,R,L,U,R) 表示处理 y=\dfrac{Px+R}{Q}\ (x\in (0,L]) 这个线段,U,R 分别为遇到横线和竖线后执行的操作序列。
既然这个东西叫万能欧几里得,自然就和欧几里得算法的迭代有关,这分为两个部分。
第一部分:P\ge Q 时,我们尝试将这个问题分解成规模为 P\bmod Q 的子问题。
我们发现,如果 P\ge Q,那么只要令 R'=U^{\frac{P}{Q}}R,然后直接递归 solve(P\bmod Q,Q,R,L,U,R') 即可,这是因为每一个 R 的前面必然有至少 \dfrac{P}{Q} 个连续的 U,它们可以和这个 R 合到一起。
第二部分:P<Q 时,交换 P,Q 地位,再进入第一部分的迭代。
容易发现,第 a 个 R 之前会有 \lfloor\dfrac{Pa+R}{Q}\rfloor 个 U,我们设第 a 个 R 在第 b 个 U 之前,那么:
b > \lfloor\dfrac{Pa+R}{Q}\rfloor
b > \dfrac{Pa+R}{Q}
a < \dfrac{Qb-R}{P}
a < \lceil\dfrac{Qb-R}{P}\rceil
因此,第 b 个 U 之前会有 \lfloor\dfrac{Qb-R-1}{P}\rfloor 个 R,这里就是一种转换思想,将 U,R 操作进行的转换。
迭代的子问题应该与原问题形式相同,所以开头的一段(第一个 U 之前)操作应该被截掉,对第一个 U 之后的部分进行递归。
假设现在去掉了第一个 U 以及它之前的 R,看看剩余的第 b 个 U 与第 b-1 个 U 之间剩余的 R 的个数:
\lfloor\dfrac{Q(b+1)-R-1}{P}\rfloor-\lfloor\dfrac{Qb-R-1}{P}\rfloor=\lfloor\dfrac{Qb+Q-R-1}{P}\rfloor-\lfloor\dfrac{Q(b-1)+Q-R-1}{P}\rfloor
这个形式就与原问题相同了,只需递归 solve(Q,P,Q-R-1,CntU-1,R,U) 即可,其中 CntU 为当前这一层 U 的总数,就等于 \lfloor\dfrac{PL+R}{Q}\rfloor,要 -1 是因为我们截去了第一个 U 之前的部分。
再考虑首位的零碎部分,首端(也就是第一个 U 之前)有 \lfloor\dfrac{Q-R-1}{P}\rfloor 个 R,所以前面要补的操作序列是 R^{\lfloor\frac{Q-R-1}{P}\rfloor}U,另外操作序列最后还有若干个 R,个数为 L-\lfloor\dfrac{Q\times CntU-R-1}{P}\rfloor,这一段也要拼到最后。
我们来看看这个过程的时间复杂度,假设一次合并操作的复杂度为 O(c),那么形如 A^{b} 这样的式子可以使用快速幂优化到 O(c\log b)。
考虑算法过程中每层迭代用到了一次 R^{\lfloor\frac{Q-R-1}{P}\rfloor},时间复杂度是 O(c \log \dfrac{Q}{P}) 的,也就是 O(c(\log Q-\log P))。
而我们又发现相邻两层的 P,Q 地位交换,所以下一层的 Q 就是这个 P,下一层的 \log Q 也就和这一层的 -\log P 抵消了。
所以总的复杂度就是 O(c\log \max(P,Q))。
Node solve (ll p,ll q,ll r,ll l,Node a,Node b) {
if (!l) {return Node();}
if (p>=q) {return solve(p%q,q,r,l,a,qpow(a,p/q)+b);}
ll m=div(l,p,r,q);
if (!m) {return qpow(b,l);}
ll cnt=l-div(q,m,-r-1,p);
return qpow(b,(q-r-1)/p)+a+solve(q,p,(q-r-1)%p,m-1,b,a)+qpow(b,cnt);
}
这题怎么做?
我们只需要讨论这题中 U 和 R 分别是什么,即可套用上面的做法。
首先考虑每个状态维护的值,我们设 s=\dfrac{Px+R}{Q},我们要求的是 s 的和,s 的平方和以及 s\times x 的和。
一般这样的题目中,还需要额外维护的是这一段操作序列中 U 和 R 的个数,此外,为了计算 s\times x 的和,我们还需要记录 x 的和。
因此我们记录六个量:U 的个数 cntu,R 的个数 cntr,x 的和 sumi,s 的和 sums,s 的平方和 sqrs,s\times x 的和 prod。
遇到 U 时,我们需要执行:cntu\leftarrow cntu+1
遇到 R 时,我们需要执行:cntr\leftarrow cntr+1,\ \ sumi\leftarrow sumi+cntr,\ \ sums\leftarrow sums+cntu,
不难推出初始状态:$U$ 的状态只有 $cntu=1$,$R$ 的状态有 $cntr,sums,sqrs$ 都为 $1$。
考虑合并状态 $S,T$,变为 $W=S+T$:
$W.cntu=S.cntu+T.cntu$,显然;
$W.cntr=S.cntr+T.cntr$,显然;
$W.sumi=S.sumi+T.sumi+S.cntr\times T.cntr$,这是因为 $T$ 中每次加的 $x$ 都比原来增大了 $S.cntr$(以左边为基数了);
$W.sums=S.sums+T.sums+S.cntu\times T.cntr$,同上;
$W.sqrs=S.sqrs+T.sqrs+T.cntr\times S.cntu^2+2\times S.cntu\times T.sums$,根据 $(s+S.cntu)^2$ 展开即可得到;
$W.prod=S.prod+T.prod+S.cntu\times S.cntr\times T.cntr+S.cntu\times T.sumi+S.cntr\times T.sums$,根据 $(s+S.cntu)(x+S.cntr)$ 展开即可得到。
最后别忘了,之前的讨论是基于 $x\in (0,L]$,本题 $x\in [0,L]$,所以最后 $sums$ 上要加一个 $\dfrac{R}{Q}$,而 $sqrs$ 上要加一个 $\dfrac{R^2}{Q^2}$。
```cpp
#include <bits/stdc++.h>
#define ll long long
using namespace std;
const int MAXN=22,P=998244353;
ll t,p,q,r,l;
struct Node {
Node () {cntu=cntr=sumi=sums=sqrs=prod=0;}
ll cntu,cntr,sumi,sums,sqrs,prod;
Node operator + (Node b) {
Node c;
c.cntu=(cntu+b.cntu)%P,c.cntr=(cntr+b.cntr)%P;
c.sumi=(sumi+b.sumi+cntr*b.cntr)%P;
c.sums=(sums+b.sums+cntu*b.cntr)%P;
c.sqrs=(sqrs+b.sqrs+((cntu*cntu)%P)*b.cntr+(2*cntu*b.sums)%P)%P;
c.prod=((prod+b.prod+((cntu*cntr)%P)*b.cntr)%P+cntu*b.sumi+cntr*b.sums)%P;
return c;
}
}nu,nr,ans;
Node qpow (Node a,ll k) {
Node res;
while (k) {
if (k&1) {res=res+a;}
a=a+a,k>>=1;
}
return res;
}
ll div (ll a,ll b,ll c,ll d) {return ((long double)1.0*a*b+c)/d;}
Node solve (ll p,ll q,ll r,ll l,Node a,Node b) {
if (!l) {return Node();}
if (p>=q) {return solve(p%q,q,r,l,a,qpow(a,p/q)+b);}
ll m=div(l,p,r,q);
if (!m) {return qpow(b,l);}
ll cnt=l-div(q,m,-r-1,p);
return qpow(b,(q-r-1)/p)+a+solve(q,p,(q-r-1)%p,m-1,b,a)+qpow(b,cnt);
}
int main () {
scanf("%lld",&t);
for (int ii=1;ii<=t;ii++) {
scanf("%lld%lld%lld%lld",&l,&p,&r,&q);
nu.cntu=1,nu.cntr=0,nu.sumi=0,nu.sums=0,nu.sqrs=0,nu.prod=0;
nr.cntu=0,nr.cntr=1,nr.sumi=1,nr.sums=0,nr.sqrs=0,nr.prod=0;
ans=qpow(nu,r/q)+solve(p,q,r%q,l,nu,nr);
printf("%lld %lld %lld\n",(ans.sums+r/q)%P,(ans.sqrs+((r/q)%P)*((r/q)%P))%P,ans.prod);
}
return 0;
}
```
------
### 这个做法好在哪里?
最关键的好处相信通过这道题目的例子大家已经了解,就是可以套板子,无论这个题长什么样,你只需要修改 $Node$ 及其加法部分就可以,而 $solve$ 部分对于所有这类题目来说都是**一模一样**的。
而且这个方法推式子没那么烦。
------
另一道万能欧几里得练手题:[LOJ #6440](https://loj.ac/problem/6440)