「算法笔记」类欧几里德算法

类欧几里德算法基于欧几里德算法,本质是通过取模缩小问题规模,递归求解。

1. $f(a,b,c,n)$

定义

推导

  1. 当 $a\ge c$ 或 $b\ge c​$ 时:

    式子 $(1)​$:我们把 $a​$ 拆成 $\left\lfloor\frac{a}{c}\right\rfloor\cdot c​$ 和 $a\bmod c​$ 两部分,$b​$ 同理。就可以将原式拆成两项了。

    式子 $(2)$:第一项直接 $O(1)$ 计算,第二项递归计算

  2. 当 $a<c​$ 且 $b<c​$ 时:

    我们将 $\frac{ai+b}{c}​$ 看作是一条线段 $y=\frac{ax+b}{c}(x\in [0,n])​$,那么原式的本质就是求线段下方、$x​$ 轴上方有多少个整点(包括线段,不包括 $x​$ 轴)。

    设 $m=\left\lfloor\frac{an+b}{c}\right\rfloor$,枚举所有点并判断是否在线段下方。

    式子 $(1)​$:转化为线段下方的整点个数。

    式子 $(2)​$:移项后再不等式右侧减去 $1​$ 使得 $\ge​$ 变为 $>​$ 号。

    式子 $(3)$:将 $>$ 一个实数转化为 $\ge$ 一个整数。我们发现如果式子 $(2)$ 中不进行转化,那么就会出现 $i\ge \frac{cj-b}{a}$,这个式子是不能直接上取整的(比如 $\lceil3\rceil=4​$ 显然是一个反例)。

    式子 $(4)​$:我们把 $\sum_{i=1}^n​$ 给去掉,直接变为一个整数。

    式子 $(5)$:将 $\sum$ 给拆开,并改变求和下界和上界

    式子 $(6)$:直接将 $\sum$ 转化为 $f$ 函数,递归计算

  3. 边界条件 $a=0$:答案为 $(n+1)\left\lfloor\frac{b}{c}\right\rfloor$。

伪代码

1
2
3
4
5
6
int f(int a,int b,int c,int n) {
if(!a) return b/c*(n+1);
if(a>=c||b>=c) return a/c*n*(n+1)/2+b/c*(n+1)+f(a%c,b%c,c,n);
int m=(a*n+b)/c;
return n*m-f(c,c-b-1,a,m-1);
}

2. $g(a,b,c,n)$

定义

推导

  1. 当 $a\ge c​$ 或 $b\ge c​$ 时:

    式子 $(1)$:我们还是仿照 $f$ 的过程拆成两项然后分别求和。

    式子 $(2)$:根据平方和公式 $\sum_{i=1}^n i^2=\frac{n(n+1)(2n+1)}{6}$,我们可以 $O(1)$ 计算出第一项的值,第二项递归计算

  2. 当 $a<c​$ 且 $b<c​$ 时:

    总体思路还是和 $f​$ 函数一样,利用其几何意义来求解,设 $m=\left\lfloor\frac{an+b}{c}\right\rfloor​$。

    式子 $(1)$:转化为线段下方的整点个数。

    式子 $(2)$:移项后再不等式右侧减去 $1$ 使得 $\ge$ 变为 $>$ 号。

    式子 $(3)$:依据还是和 $f$ 函数相同。

    式子 $(4)​$:由于第二个 $\sum​$ 内就是对 $i\in \left[\left\lfloor\frac{cj-b-1}{a}\right\rfloor+1,n\right]​$ 求和,因此直接利用等差数列求和

    式子 $(5)$:将等差数列求和公式继续展开,转化为 $3​$ 项。

    式子 $(6)​$:第一项可以 $O(1)​$ 计算,第二项和第三项发现是 $f​$ 和 $h​$ 函数,可以递归计算

  3. 边界条件 $a=0​$:答案为 $\frac{n(n+1)}{2}\left\lfloor\frac{b}{c}\right\rfloor​$。

伪代码

1
2
3
4
5
6
int g(int a,int b,int c,int n) {
if(!a) return b/c*n*(n+1)/2;
if(a>=c||b>=c) return a/c*n*(n+1)*(2*n+1)/6+b/c*n*(n+1)/2+g(a%c,b%c,c,n);
int m=(a*n+b)/c;
return m*n*(n+1)/2-f(c,c-b-1,a,m-1)/2-h(c,c-b-1,a,m-1)/2;
}

3. $h(a,b,c,n)$

定义

推导

  1. 当 $a\ge c​$ 或 $b\ge c​$ 时:

    式子 $(1)$:把原式的 $a$ 和 $b$ 均拆成两个部分。

    式子 $(2)​$:暴力展开平方项。

    式子 $(3)$:将展开式中能 $O(1)$ 直接计算的提出,剩下的转化为 $f,g,h​$ 函数递归计算

  2. 当 $a<c$ 且 $b<c$ 时:

    我们试图把平方项用另一种形式表示出来:$n^2=2\sum_{i=1}^n i-n$。我们就根据这个公式进行转化(其中 $m$ 的值依旧是 $\left\lfloor\frac{an+b}{c}\right\rfloor​$。

    式子 $(1)$:我们按照上述公式将 $\left\lfloor\frac{ai+b}{c}\right\rfloor^2$ 进行变形。

    式子 $(2)$:变换求和顺序,将 $j$ 的上界限定在 $\left\lfloor\frac{ai+b}{c}\right\rfloor$ 即 $\frac{ai+b}{c}$ 范围内。上述公式中多出来的一项转化为 $f​$ 函数。

    式子 $(3)$:改变 $j$ 的上下界,并按照前文的方法(详见 $f$ 或 $g$ 函数的这部分内容)确定 $i$ 的范围。

    式子 $(4)​$:将求和式展开。

    式子 $(5)​$:将能够直接计算的提出,其余递归计算

  3. 边界条件 $a=0$:答案为 $(n+1)\left\lfloor\frac{b}{c}\right\rfloor^2$。

伪代码

1
2
3
4
5
6
int h(int a,int b,int c,int n) {
if(!a) return sqr(b/c)*(n+1);
if(a>=c||b>=c) return sqr(a/c)*n*(n+1)*(2*n+1)/6+sqr(b/c)*(n+1)+(a/c)*(b/c)*n*(n+1)+h(a%c,b%c,c,n)+2*(a/c)*g(a%c,b%c,c,n)+2*(b/c)*f(a%c,b%c,c,n);
int m=(a*n+b)/c;
return n*m*(m+1)-f(a,b,c,n)-2*g(c,c-b-1,a,m-1)-2*f(c,c-b-1,a,m-1);
}

整体求解

如果我们要对 $f,g,h$ 三个函数同时求解,显然是没法记忆化的。

注意到每个函数中都调用了 $(a\bmod c,b\bmod c,c,n)$ 和 $(c,c-b-1,a,m-1)$,那么我们考虑将 $3$ 个函数在同一个函数中同时求解

时间复杂度:$O(\log n)$

代码如下:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
#include <cstdio>
typedef long long LL;

const int mod=998244353;
const LL i2=499122177,i6=166374059;

struct Ans {
LL f,g,h;
};
LL sqr(LL x) {
return x*x%mod;
}
LL S2(LL n) {
return n*(n+1)%mod*i2%mod;
}
LL S3(LL n) {
return n*(n+1)%mod*(n+n+1)%mod*i6%mod;
}
Ans solve(LL a,LL b,LL c,LL n) {
Ans now;
LL A=a/c,B=b/c;
if(!a) {
now.f=(n+1)*B%mod;
now.g=S2(n)*B%mod;
now.h=(n+1)*sqr(B)%mod;
return now;
}
if(a>=c||b>=c) {
Ans nxt=solve(a%c,b%c,c,n);
now.f=S2(n)*A%mod+(n+1)*B%mod+nxt.f;
now.g=S3(n)*A%mod+S2(n)*B%mod+nxt.g;
now.h=S3(n)*sqr(A)%mod+(n+1)*sqr(B)%mod+2*S2(n)*A%mod*B%mod+(B+B)*nxt.f%mod+(A+A)*nxt.g%mod+nxt.h;
now.f%=mod,now.g%=mod,now.h%=mod;
return now;
}
LL m=(a*n+b)/c;
Ans nxt=solve(c,c-b-1,a,m-1);
now.f=n*m%mod-nxt.f;
now.g=(m*n%mod*(n+1)%mod-nxt.f-nxt.h)*i2%mod;
now.h=n*m%mod*(m+1)%mod-now.f-2*(nxt.f+nxt.g);
now.f%=mod,now.g%=mod,now.h%=mod;
return now;
}
int main() {
int T;
for(scanf("%d",&T);T--;) {
LL a,b,c,n;
scanf("%lld%lld%lld%lld",&n,&a,&b,&c);
Ans ans=solve(a,b,c,n);
ans.f=(ans.f%mod+mod)%mod;
ans.g=(ans.g%mod+mod)%mod;
ans.h=(ans.h%mod+mod)%mod;
printf("%lld %lld %lld\n",ans.f,ans.g,ans.h);
}
return 0;
}

习题


附:高次方和求法

首先,我们有如下公式:

接下来我们就分析一下这些公式是怎么得到的。

根据

我们可以得到

具体证明只需要将右侧式子按照组合数的基本公式 $\binom{n}{k}=\binom{n-1}{k-1}+\binom{n-1}{k}​$ 暴力展开。

那么我们可以根据这个结论得到:

对于平方求和,证明如下:

再证明一下立方求和:

更高的次数也能按照以上规则写出求和式!

如果用代码实现的话,可以采用拉格朗日插值(详见「算法笔记」拉格朗日插值)求出 $k$ 次方和,时间复杂度最优可以达到 $O(k)$。

0%