【LOJ138】类欧几里得算法

最新推荐文章于 2023-11-21 12:28:31 发布

cz_xuyixuan

最新推荐文章于 2023-11-21 12:28:31 发布

阅读量947

点赞数 3

分类专栏：【类型】做题记录【OJ】LOJ 【算法】拉格朗日插值法【算法】欧几里得算法【算法】类欧几里得算法【资料】模板题

本文链接：https://blog.csdn.net/qq_39972971/article/details/94618394

版权

【类型】做题记录同时被 3 个专栏收录

767 篇文章 2 订阅

订阅专栏

【OJ】LOJ

123 篇文章 0 订阅

订阅专栏

【资料】模板题

54 篇文章 0 订阅

订阅专栏

【题目链接】

点击打开链接

【思路要点】

以下考虑实现函数 $f u n c (N, a, b, c)$ ，计算 $0\leq k_1+k_2\leq10$ 的情况下所求式子的值，即
$\sum_{i=0}^{N}i^{k_1}\lfloor\frac{ai+b}{c}\rfloor^{k_2}$
若 $a = 0$ 或 $\lfloor\frac{aN+b}{c}\rfloor=0$ ，那么， $\lfloor\frac{ai+b}{c}\rfloor$ 的取值始终相同，答案即为
$\lfloor\frac{aN+b}{c}\rfloor\sum_{i=0}^{N}i^{k_1}$
可以用插值解决。
若 $a\geq c$ ，则令 $q=\lfloor\frac{a}{c}\rfloor,r=a-qc$ ，所求式子即为
$\sum_{i=0}^{N}i^{k_1}(qi+\lfloor\frac{ri+b}{c}\rfloor)^{k_2}\\=\sum_{i=0}^{N}i^{k_1}\sum_{j=0}^{k_2}\binom{k_2}{j}(qi)^j\lfloor\frac{ri+b}{c}\rfloor^{k_2-j}\\=\sum_{j=0}^{k_2}\binom{k_2}{j}q^j\sum_{i=0}^{N}i^{k_1+j}\lfloor\frac{ri+b}{c}\rfloor^{k_2-j}$
递归计算 $f u n c (N, r, b, c)$ 即可。
若 $b\geq c$ ，则令 $q=\lfloor\frac{b}{c}\rfloor,r=b-qc$ ，所求式子即为
$\sum_{i=0}^{N}i^{k_1}(q+\lfloor\frac{ai+r}{c}\rfloor)^{k_2}\\=\sum_{i=0}^{N}i^{k_1}\sum_{j=0}^{k_2}\binom{k_2}{j}q^j\lfloor\frac{ai+r}{c}\rfloor^{k_2-j}\\=\sum_{j=0}^{k_2}\binom{k_2}{j}q^j\sum_{i=0}^{N}i^{k_1}\lfloor\frac{ai+r}{c}\rfloor^{k_2-j}$
递归计算 $f u n c (N, a, r, c)$ 即可。
否则，即 $a<c,b<c,\lfloor\frac{aN+b}{c}\rfloor>0$ ，令 $M=\lfloor\frac{aN+b}{c}\rfloor$ 。
注意到 $\lfloor\frac{ai+b}{c}\rfloor^{k_2}$ 较难处理，需要将其变形，有
$\lfloor\frac{ai+b}{c}\rfloor^{k_2}=\sum_{j=0}^{\lfloor\frac{ai+b}{c}\rfloor-1}((j+1)^{k_2}-j^{k_2})$
注意这里规定 $0^0=0$ 。
上述变换的好处是提供了交换求和顺序的可能，由此，原式可变为
$\sum_{i=0}^{N}i^{k_1}\sum_{j=0}^{\lfloor\frac{ai+b}{c}\rfloor-1}((j+1)^{k_2}-j^{k_2})\\=\sum_{j=0}^{M-1}((j+1)^{k_2}-j^{k_2})\sum_{i=0}^{N}i^{k_1}[\lfloor\frac{ai+b}{c}\rfloor\geq j+1]\\=\sum_{j=0}^{M-1}((j+1)^{k_2}-j^{k_2})\sum_{i=0}^{N}i^{k_1}[i>\lfloor\frac{cj+c-b-1}{a}\rfloor]\\=\sum_{j=0}^{M-1}((j+1)^{k_2}-j^{k_2})\sum_{i=0}^{N}i^{k_1}-\sum_{j=0}^{M-1}((j+1)^{k_2}-j^{k_2})\sum_{i=0}^{\lfloor\frac{cj+c-b-1}{a}\rfloor}i^{k_1}$
靠前的求和符号可以直接计算，考虑靠后的求和符号。
$j+1)^{k_2}-j^{k_2}$ 显然是关于 $j$ 的 $k_2-1$ 次多项式，而 $\sum_{i=0}^{\lfloor\frac{cj+c-b-1}{a}\rfloor}i^{k_1}$ 也是关于 $\lfloor\frac{cj+c-b-1}{a}\rfloor$ 的 $k_1+1$ 次多项式，不妨令为 $A (x), B (x)$ 。
那么
$\sum_{j=0}^{M-1}((j+1)^{k_2}-j^{k_2})\sum_{i=0}^{\lfloor\frac{cj+c-b-1}{a}\rfloor}i^{k_1}\\=\sum_{i=0}^{M-1}\sum_{j=0}^{k_2-1}A_ii^j\sum_{k=0}^{k_1+1}B_k\lfloor\frac{cj+c-b-1}{a}\rfloor^{k}\\=\sum_{j=0}^{k_2-1}A_i\sum_{k=0}^{k_1+1}B_k\sum_{i=0}^{M-1}i^j\lfloor\frac{cj+c-b-1}{a}\rfloor^{k}$
递归计算 $f u n c (M - 1, c, c - b - 1, a)$ 即可。
注意到在递归中 $(a,c)\rightarrow(a\%c,c)\rightarrow(c,a\%c)$ ，因此递归层数为 $O (L o g V)$ 。
时间复杂度 $O(TK^4LogV)$ 。

【代码】

#include<bits/stdc++.h>
using namespace std;
const int MAXN = 15;
const int P = 1e9 + 7;
typedef long long ll;
typedef long double ld;
typedef unsigned long long ull;
template <typename T> void chkmax(T &x, T y) {x = max(x, y); }
template <typename T> void chkmin(T &x, T y) {x = min(x, y); } 
template <typename T> void read(T &x) {
	x = 0; int f = 1;
	char c = getchar();
	for (; !isdigit(c); c = getchar()) if (c == '-') f = -f;
	for (; isdigit(c); c = getchar()) x = x * 10 + c - '0';
	x *= f;
}
template <typename T> void write(T x) {
	if (x < 0) x = -x, putchar('-');
	if (x > 9) write(x / 10);
	putchar(x % 10 + '0');
}
template <typename T> void writeln(T x) {
	write(x);
	puts("");
}
struct info {
	int a[MAXN][MAXN];
};
int sum[MAXN][MAXN];
int binom[MAXN][MAXN];
int power(int x, int y) {
	if (y == 0) return 1;
	int tmp = power(x, y / 2);
	if (y % 2 == 0) return 1ll * tmp * tmp % P;
	else return 1ll * tmp * tmp % P * x % P;
}
void update(int &x, int y) {
	x += y;
	if (x >= P) x -= P;
}
info func(int n, int a, int b, int c) {
	assert(n >= 0 && a >= 0 && b >= 0 && c >= 0);
	info ans;
	memset(ans.a, 0, sizeof(ans.a));
	if (a == 0 || (1ll * a * n + b) / c == 0) {
		for (int k1 = 0; k1 <= 10; k1++) {
			int mul = 0, now = 1;
			for (int i = 0; i <= k1 + 1; i++) {
				update(mul, 1ll * now * sum[k1][i] % P);
				now = 1ll * now * n % P;
			}
			int base = (1ll * a * n + b) / c % P; now = 1;
			for (int k2 = 0; k1 + k2 <= 10; k2++) {
				ans.a[k1][k2] = 1ll * now * mul % P;
				now = 1ll * now * base % P;
			}
		}
		return ans;
	}
	if (a >= c) {
		info tmp = func(n, a % c, b, c);
		for (int k1 = 0; k1 <= 10; k1++)
		for (int k2 = 0; k1 + k2 <= 10; k2++) {
			int now = 1, base = a / c;
			for (int i = 0; i <= k2; i++) {
				update(ans.a[k1][k2], 1ll * binom[k2][i] * now % P * tmp.a[k1 + i][k2 - i] % P);
				now = 1ll * now * base % P;
			}
		}
		return ans;
	}
	if (b >= c) {
		info tmp = func(n, a, b % c, c);
		for (int k1 = 0; k1 <= 10; k1++)
		for (int k2 = 0; k1 + k2 <= 10; k2++) {
			int now = 1, base = b / c;
			for (int i = 0; i <= k2; i++) {
				update(ans.a[k1][k2], 1ll * binom[k2][i] * now % P * tmp.a[k1][k2 - i] % P);
				now = 1ll * now * base % P;
			}
		}
		return ans;
	}
	int m = (1ll * a * n + b) / c;
	info tmp = func(m - 1, c, c - b - 1, a);
	for (int k1 = 0; k1 <= 10; k1++) {
		int all = 0, now = 1;
		for (int i = 0; i <= k1 + 1; i++) {
			update(all, 1ll * now * sum[k1][i] % P);
			now = 1ll * now * n % P;
		}
		for (int k2 = 0; k1 + k2 <= 10; k2++) {
			ans.a[k1][k2] = 1ll * power(m, k2) * all % P;
			for (int i = 0; i <= k2 - 1; i++)
			for (int j = 0; j <= k1 + 1; j++) {
				int coef = 1ll * binom[k2][i] * sum[k1][j] % P;
				update(ans.a[k1][k2], P - 1ll * coef * tmp.a[i][j] % P);
			}
		}
	}
	return ans;
}
void Lagrange(int n, int *x, int *y, int *a) {
	static int p[MAXN], q[MAXN];
	memset(p, 0, sizeof(p)); p[0] = 1;
	for (int i = 1; i <= n; i++) {
		for (int j = i - 1; j >= 0; j--) {
			p[j + 1] = (p[j + 1] + p[j]) % P;
			p[j] = (P - 1ll * p[j] * x[i] % P) % P;
		}
	}
	for (int i = 1; i <= n; i++) {
		memset(q, 0, sizeof(q));
		for (int j = n - 1; j >= 0; j--)
			q[j] = (p[j + 1] + 1ll * q[j + 1] * x[i]) % P;
		int now = 1;
		for (int j = 1; j <= n; j++)
			if (j != i) now = 1ll * now * (x[i] - x[j]) % P;
		now = power((P + now) % P, P - 2);
		for (int j = 0; j <= n; j++)
			q[j] = 1ll * q[j] * now % P;
		for (int j = 0; j <= n; j++)
			a[j] = (a[j] + 1ll * q[j] * y[i]) % P;
	}
}
int main() {
	for (int i = 0; i <= 10; i++) {
		static int pos[MAXN], now[MAXN];
		now[0] = power(0, i);
		for (int j = 1; j <= i + 2; j++) {
			now[j] = (now[j - 1] + power(j, i)) % P;
			pos[j] = j;
		}
		Lagrange(i + 2, pos, now, sum[i]);
	}
	for (int i = 0; i <= 10; i++) {
		binom[i][0] = 1;
		for (int j = 1; j <= i; j++)
			binom[i][j] = binom[i - 1][j - 1] + binom[i - 1][j];
	}
	int T; read(T);
	while (T--) {
		int n, a, b, c, k1, k2;
		read(n), read(a), read(b), read(c), read(k1), read(k2);
		writeln(func(n, a, b, c).a[k1][k2]);
	}
	return 0;
}