为了简化代码,定义以下函数,下文不再重复:
$$ \begin{aligned} S0(n) &= \sum_{i = 0} ^ n i ^ 0 = n + 1 \\ S1(n) &= \sum_{i = 0} ^ n i ^ 1 = \frac{n(n + 1)}{2} \\ S2(n) &= \sum_{i = 0} ^ n i ^ 2 = \frac{n(n+1)(2n+1)}{6} \end{aligned} $$
1. $f(a,b,c,n)$
定义
$$ f(a,b,c,n)=\sum_{i=0}^n \left\lfloor\frac{ai+b}{c}\right\rfloor $$
推导
当 $a\ge c$ 或 $b\ge c$ 时:
$$ \begin{align*} f(a,b,c,n) & = \sum_{i=0}^n \left\lfloor\frac{ai+b}{c}\right\rfloor \\ & = \sum_{i=0}^n \left(i\left\lfloor\frac{a}{c}\right\rfloor+\left\lfloor\frac{b}{c}\right\rfloor\right)+\sum_{i=0}^n \left\lfloor\frac{a\bmod c\times i+b\bmod c}{c}\right\rfloor\tag 1 \\ & = \frac{n(n+1)}{2}\left\lfloor\frac{a}{c}\right\rfloor+(n+1)\left\lfloor\frac{b}{c}\right\rfloor + f(a\bmod c,b\bmod c,c,n)\tag 2 \end{align*} $$
式子 $(1)$:我们把 $a$ 拆成 $\left\lfloor\frac{a}{c}\right\rfloor\cdot c$ 和 $a\bmod c$ 两部分,$b$ 同理。就可以将原式拆成两项了。
式子 $(2)$:第一项直接 $\mathcal O(1)$ 计算,第二项递归计算。
当 $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$,枚举所有点并判断是否在线段下方。
$$ \begin{align*} f(a,b,c,n) & = \sum_{i=0}^n \left\lfloor\frac{ai+b}{c}\right\rfloor \\ & = \sum_{i=0}^n \sum_{j=1}^m \left[j\le\frac{ai+b}{c}\right]\tag 1 \\ & = \sum_{i=0}^n \sum_{j=1}^m \left[ai > cj-b-1\right]\tag 2 \\ & = \sum_{j=1}^m \sum_{i=0}^n \left[i\ge \left\lfloor\frac{cj-b-1}{a}\right\rfloor+1\right]\tag 3 \\ & = \sum_{j=1}^m \left(n-\left\lfloor\frac{cj-b-1}{a}\right\rfloor\right)\tag 4 \\ & = nm-\sum_{j=0}^{m-1} \left\lfloor\frac{cj+c-b-1}{a}\right\rfloor\tag 5 \\ & = nm-f(c,c-b-1,a,m-1)\tag 6 \end{align*} $$
式子 $(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$ 函数,递归计算。
- 边界条件 $a=0$:答案为 $(n+1)\left\lfloor\frac{b}{c}\right\rfloor$。
伪代码
int f(int a, int b, int c, int n) {
if (!a) {
return b / c * S0(n);
}
if (a >= c || b >= c) {
return a / c * S1(n) + b / c * S0(n) + 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)$
定义
$$ g(a,b,c,n)=\sum_{i=0}^n i\left\lfloor\frac{ai+b}{c}\right\rfloor $$
推导
当 $a\ge c$ 或 $b\ge c$ 时:
$$ \begin{align*} g(a,b,c,n) & =\sum_{i=0}^n i\left\lfloor\frac{ai+b}{c}\right\rfloor \\ & =\sum_{i=0}^n \left(i^2\left\lfloor\frac{a}{c}\right\rfloor+i\left\lfloor\frac{b}{c}\right\rfloor\right)+\sum_{i=0}^n i\left\lfloor\frac{a\bmod c\times i+b\bmod c}{c}\right\rfloor\tag 1 \\ & =\frac{n(n+1)(2n+1)}{6}\left\lfloor\frac{a}{c}\right\rfloor+\frac{n(n+1)}{2}\left\lfloor\frac{b}{c}\right\rfloor+g(a\bmod c,b\bmod c,c,n)\tag 2 \end{align*} $$
式子 $(1)$:我们还是仿照 $f$ 的过程拆成两项然后分别求和。
式子 $(2)$:根据平方和公式 $\sum_{i=1}^n i^2=\frac{n(n+1)(2n+1)}{6}$,我们可以 $\mathcal O(1)$ 计算出第一项的值,第二项递归计算。
当 $a < c$ 且 $b < c$ 时:
总体思路还是和 $f$ 函数一样,利用其几何意义来求解,设 $m=\left\lfloor\frac{an+b}{c}\right\rfloor$。
$$ \begin{align*} g(a,b,c,n) & = \sum_{i=0}^n i\left\lfloor\frac{ai+b}{c}\right\rfloor \\ & = \sum_{i=0}^n i\times \sum_{j=1}^m \left[j\le\frac{ai+b}{c}\right]\tag 1 \\ & = \sum_{i=0}^n i\times \sum_{j=1}^m \left[ai > cj-b-1\right]\tag 2 \\ & = \sum_{j=1}^m \sum_{i=0}^n i\left[i\ge \left\lfloor\frac{cj-b-1}{a}\right\rfloor+1\right]\tag 3 \\ & = \sum_{j=1}^m \frac{1}{2}\left(n+\left\lfloor\frac{cj-b-1}{a}\right\rfloor+1\right)\left(n-\left\lfloor\frac{cj-b-1}{a}\right\rfloor\right)\tag 4 \\ & = \frac{mn(n+1)}{2}-\frac{1}{2}\sum_{j=0}^{m-1}\left\lfloor\frac{cj+c-b-1}{a}\right\rfloor-\frac{1}{2}\sum_{j=0}^{m-1}\left\lfloor\frac{cj+c-b-1}{a}\right\rfloor^2 \tag 5 \\ & = \frac{mn(n+1)}{2}-\frac{1}{2}f(c,c-b-1,a,m-1)-\frac{1}{2}h(c,c-b-1,a,m-1) \tag 6 \end{align*} $$
式子 $(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)$:第一项可以 $\mathcal O(1)$ 计算,第二项和第三项发现是 $f$ 和 $h$ 函数,可以递归计算。
- 边界条件 $a=0$:答案为 $\frac{n(n+1)}{2}\left\lfloor\frac{b}{c}\right\rfloor$。
伪代码
int g(int a, int b, int c, int n) {
if (!a) {
return b / c * S1(n);
}
if (a >= c || b >= c) {
return a / c * S2(n) + b / c * S1(n) + g(a % c, b % c, c, n);
}
int m = (a * n + b) / c;
return m * S1(n) - f(c, c - b - 1, a, m - 1) / 2 - h(c, c - b - 1, a, m - 1) / 2;
}
3. $h(a,b,c,n)$
定义
$$ h(a,b,c,n)=\sum_{i=0}^n \left\lfloor\frac{ai+b}{c}\right\rfloor^2 $$
推导
当 $a\ge c$ 或 $b\ge c$ 时:
$$ \begin{align*} h(a,b,c,n) & =\sum_{i=0}^n \left\lfloor\frac{ai+b}{c}\right\rfloor^2 \\ & = \sum_{i=0}^n \left(i\left\lfloor\frac{a}{c}\right\rfloor+\left\lfloor\frac{b}{c}\right\rfloor+\left\lfloor\frac{a\bmod c\times i+b\bmod c}{c}\right\rfloor\right)^2 \tag 1 \\ & = \sum_{i=0}^n i^2\left\lfloor\frac{a}{c}\right\rfloor^2+\left\lfloor\frac{b}{c}\right\rfloor^2+\left\lfloor\frac{a\bmod c\times i+b\bmod c}{c}\right\rfloor^2+2i\left\lfloor\frac{a}{c}\right\rfloor\left\lfloor\frac{b}{c}\right\rfloor+2i\left\lfloor\frac{a}{c}\right\rfloor\left\lfloor\frac{a\bmod c\times i+b\bmod c}{c}\right\rfloor+2\left\lfloor\frac{b}{c}\right\rfloor\left\lfloor\frac{a\bmod c\times i+b\bmod c}{c}\right\rfloor \tag 2\\ & = \frac{n(n+1)(2n+1)}{6}\left\lfloor\frac{a}{c}\right\rfloor^2+(n+1)\left\lfloor\frac{b}{c}\right\rfloor^2+n(n+1)\left\lfloor\frac{a}{c}\right\rfloor\left\lfloor\frac{b}{c}\right\rfloor+h(a\bmod c,b\bmod c,c,n)+2\left\lfloor\frac{a}{c}\right\rfloor g(a\bmod c,b\bmod c,c,n)+2\left\lfloor\frac{b}{c}\right\rfloor f(a\bmod c,b\bmod c,c,n) \tag 3 \end{align*} $$
式子 $(1)$:把原式的 $a$ 和 $b$ 均拆成两个部分。
式子 $(2)$:暴力展开平方项。
式子 $(3)$:将展开式中能 $\mathcal O(1)$ 直接计算的提出,剩下的转化为 $f,g,h$ 函数递归计算。
当 $a < c$ 且 $b < c$ 时:
我们试图把平方项用另一种形式表示出来:$n^2=2\sum_{i=1}^n i-n$。我们就根据这个公式进行转化(其中 $m$ 的值依旧是 $\left\lfloor\frac{an+b}{c}\right\rfloor$。
$$ \begin{align*} h(a,b,c,n) & =\sum_{i=0}^n \left\lfloor\frac{ai+b}{c}\right\rfloor^2 \\ & = \sum_{i=0}^n \left(2\times\sum_{j=1}^{\left\lfloor\frac{ai+b}{c}\right\rfloor} j-\left\lfloor\frac{ai+b}{c}\right\rfloor\right) \tag 1 \\ & = 2\times \sum_{j=1}^m \sum_{i=0}^n j\times\left[j\le \frac{ai+b}{c}\right]-f(a,b,c,n) \tag 2 \\ & = 2\times \sum_{j=0}^{m-1} (j+1)\left(n-\left\lfloor\frac{cj+c-b-1}{a}\right\rfloor\right)-f(a,b,c,n)\tag 3 \\ & = 2\times\sum_{j=0}^{m-1} \left(n(j+1)-j\left\lfloor\frac{cj+c-b-1}{a}\right\rfloor-\left\lfloor\frac{cj+c-b-1}{a}\right\rfloor\right)-f(a,b,c,n) \tag 4 \\ & = nm(m+1)-2g(c,c-b-1,a,m-1)-2f(c,c-b-1,a,m-1)-f(a,b,c,n)\tag 5 \end{align*} $$
式子 $(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)$:将能够直接计算的提出,其余递归计算。
- 边界条件 $a=0$:答案为 $(n+1)\left\lfloor\frac{b}{c}\right\rfloor^2$。
伪代码
int h(int a, int b, int c, int n) {
if (!a) {
return sqr(b / c) * S0(n);
}
if (a >= c || b >= c) {
return sqr(a / c) * S2(n) + sqr(b / c) * S0(n) + (a / c) * (b / c) * 2 * S1(n)
+ 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 * S0(m) - 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$ 个函数在同一个函数中同时求解。
时间复杂度:$\mathcal O(\log n)$
代码如下:
#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 S0(LL n) {
return n + 1;
}
LL S1(LL n) {
return n * (n + 1) % mod * i2 % mod;
}
LL S2(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 = S0(n) * B % mod;
now.g = S1(n) * B % mod;
now.h = S0(n) * sqr(B) % mod;
return now;
}
if (a >= c || b >= c) {
Ans nxt = solve(a % c, b % c, c, n);
now.f = S1(n) * A % mod + S0(n) * B % mod + nxt.f;
now.g = S2(n) * A % mod + S1(n) * B % mod + nxt.g;
now.h = S2(n) * sqr(A) % mod + S0(n) * sqr(B) % mod + 2 * S1(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 * S0(n) % mod - nxt.f - nxt.h) * i2 % mod;
now.h = n * m % mod * S0(m) % mod - now.f - 2 * (nxt.f + nxt.g);
now.f %= mod, now.g %= mod, now.h %= mod;
return now;
}
int main() {
int m;
scanf("%d", &m);
for (int i = 1; i <= m; i++) {
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;
}
高次方和求法
首先,我们有如下公式:
$$ \begin{aligned} &\sum_{i=1}^n i^1=\frac{n(n+1)}{2}\\ &\sum_{i=1}^n i^2=\frac{n(n+1)(2n+1)}{6}\\ &\sum_{i=1}^n i^3=\frac{n^2(n+1)^2}{4}\\ &\sum_{i=1}^n i^4=\cdots \end{aligned} $$
接下来我们就分析一下这些公式是怎么得到的。
根据
$$ \binom{n}{k}=\binom{n-1}{k-1}+\binom{n-1}{k} $$
我们可以得到
$$ \sum_{i=1}^n \binom{i}{k}=\binom{n+1}{k+1} $$
具体证明只需要将右侧式子按照组合数的基本公式 $\binom{n}{k}=\binom{n-1}{k-1}+\binom{n-1}{k}$ 暴力展开。
那么我们可以根据这个结论得到:
$$ \begin{align*} \sum_{i=1}^n i^1 & = \sum_{i=1}^n i=\sum_{i=1}^n \binom{i}{1}=\binom{n+1}{2} \\ & =\frac{n(n+1)}{2} \end{align*} $$
对于平方求和,证明如下:
$$ \begin{align*} \sum_{i=1}^n i^2 &= \sum_{i=1}^n \left[i(i-1)+i\right]=2\sum_{i=1}^n \binom{i}{2}+\sum_{i=1}^n \binom{i}{1}=2\binom{n+1}{3}+\binom{n+1}{2} \\ &=\frac{n(n+1)(2n+1)}{6} \end{align*} $$
再证明一下立方求和:
$$ \begin{align*} \sum_{i=1}^n i^3&=\sum_{i=1}^n \left[i(i-1)(i-2)+3i(i-1)+i\right]=6\sum_{i=1}^n\binom{i}{3}+6\sum_{i=1}^n\binom{i}{2}+\sum_{i=1}^n\binom{i}{1}\\ &=6\binom{n+1}{4}+6\binom{n+1}{3}+\binom{n+1}{2} \\ &=\frac{n^2(n+1)^2}{4} \end{align*} $$
更高的次数也能按照以上规则写出求和式!
如果用代码实现的话,可以采用拉格朗日插值(详见「算法笔记」拉格朗日插值)求出 $k$ 次方和,时间复杂度最优可以达到 $\mathcal O(\min(k, n, P))$,其中 $n$ 是数字个数,$P$ 是模数。