【LOJ3058】「HNOI2019」白兔之舞

最新推荐文章于 2019-12-28 17:15:46 发布
最新推荐文章于 2019-12-28 17:15:46 发布
阅读量592 收藏
点赞数
CC 4.0 BY-SA版权
分类专栏： 【类型】做题记录 【OJ】LOJ 【算法】FFT与NTT 【算法】单位根反演 【算法】矩阵乘法
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.youkuaiyun.com/qq_39972971/article/details/89512679
本文介绍了一种解决特定矩阵幂运算问题的方法，利用任意模数FFT进行快速计算。通过求解矩阵的kkk次单位根，结合矩阵快速幂技巧，实现了O(L+N^3KLogL+KLogK)的时间复杂度。代码示例展示了如何初始化FFT，执行前缀和变换，以及如何通过多项式乘法加速计算。

【题目链接】

【思路要点】

  • 首先，求出 L L L 的任意原根 g g g ，以及其 k k k 次单位根 w w w 。
  • 记初始给定的矩阵为 m a t mat mat ，那么答案 a n s t ans_t anst​ 应当满足
    a n s t = ∑ i = 0 L ( L i ) ( m a t i ) x , y [ i % k = t ] ans_t=\sum_{i=0}^{L}\binom{L}{i}(mat^i)_{x,y}[i\%k=t] anst​=i=0∑L​(iL​)(mati)x,y​[i%k=t]
  • 将 [ i % k = t ] [i\%k=t] [i%k=t] 展开为 1 k ∑ j = 0 k − 1 w j ( i − t ) \frac{1}{k}\sum_{j=0}^{k-1}w^{j(i-t)} k1​∑j=0k−1​wj(i−t) ，则有
    k × a n s t = ∑ i = 0 L ( L i ) ( m a t i ) x , y ∑ j = 0 k − 1 w j ( i − t ) k\times ans_t=\sum_{i=0}^{L}\binom{L}{i}(mat^i)_{x,y}\sum_{j=0}^{k-1}w^{j(i-t)} k×anst​=i=0∑L​(iL​)(mati)x,y​j=0∑k−1​wj(i−t)
    k × a n s t = ∑ j = 0 k − 1 w − j t ∑ i = 0 L ( L i ) w i j ( m a t i ) x , y k\times ans_t=\sum_{j=0}^{k-1}w^{-jt}\sum_{i=0}^{L}\binom{L}{i}w^{ij}(mat^i)_{x,y} k×anst​=j=0∑k−1​w−jti=0∑L​(iL​)wij(mati)x,y​
    k × a n s t = ∑ j = 0 k − 1 w − j t ( ( w j m a t + I ) L ) x , y k\times ans_t=\sum_{j=0}^{k-1}w^{-jt}((w^jmat+I)^L)_{x,y} k×anst​=j=0∑k−1​w−jt((wjmat+I)L)x,y​
  • 那么 w − ( t 2 ) × k × a n s t = ∑ j = 0 k − 1 w − ( j + t 2 ) × w ( j 2 ) ( ( w j m a t + I ) L ) x , y w^{-\binom{t}{2}}\times k\times ans_t=\sum_{j=0}^{k-1}w^{-\binom{j+t}{2}}\times w^{\binom{j}{2}}((w^jmat+I)^L)_{x,y} w−(2t​)×k×anst​=j=0∑k−1​w−(2j+t​)×w(2j​)((wjmat+I)L)x,y​
  • 使用任意模数 F F T FFT FFT 计算之即可。
  • 时间复杂度 O ( L + N 3 K L o g L + K L o g K ) O(\sqrt{L}+N^3KLogL+KLogK) O(L ​+N3KLogL+KLogK) 。

【代码】

#include<bits/stdc++.h>
using namespace std;
const int MAXN = 262144;
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("");
}
namespace AnyModuloFFT {
	const int MAXN = 262144;
	const long double pi = acosl(-1);
	struct point {long double x, y; };
	point operator + (point a, point b) {return (point) {a.x + b.x, a.y + b.y}; }
	point operator - (point a, point b) {return (point) {a.x - b.x, a.y - b.y}; }
	point operator * (point a, point b) {return (point) {a.x * b.x - a.y * b.y, a.x * b.y + a.y * b.x}; }
	point operator / (point a, long double x) {return (point) {a.x / x, a.y / x}; }
	int N, Log, home[MAXN];
	point tmp[MAXN];
	void FFTinit() {
		for (int i = 0; i < N; i++) {
			int tmp = i, ans = 0;
			for (int j = 1; j <= Log; j++) {
				ans <<= 1;
				ans += tmp & 1;
				tmp >>= 1;
			}
			home[i] = ans;
		}
	}
	void FFT(point *a, int mode) {
		for (int i = 0; i < N; i++)
			if (home[i] < i) swap(a[i], a[home[i]]);
		for (int len = 2; len <= N; len <<= 1) {
			point delta = (point) {cosl(2 * pi / len * mode), sinl(2 * pi / len * mode)};
			for (int i = 0; i < N; i += len) {
				point now = (point) {1, 0};
				for (int j = i, k = i + len / 2; k < i + len; j++, k++) {
					point tmp = a[j];
					point tnp = a[k] * now;
					a[j] = tmp + tnp;
					a[k] = tmp - tnp;
					now = now * delta;
				}
			}
		}
		if (mode == -1) {
			for (int i = 0; i < N; i++)
				a[i] = a[i] / N;
		}
	}
	void times(int *a, int *b, int *c, int P, int limit) {
		N = 1, Log = 0;
		while (N <= 2 * limit) {
			N <<= 1;
			Log++;
		}
		static point ax[MAXN], ay[MAXN];
		static point bx[MAXN], by[MAXN];
		for (int i = 0; i <= limit; i++) {
			ax[i] = (point) {a[i] & 32767, 0};
			ay[i] = (point) {a[i] >> 15, 0};
			bx[i] = (point) {b[i] & 32767, 0};
			by[i] = (point) {b[i] >> 15, 0};
		}
		for (int i = limit + 1; i < N; i++) {
			ax[i] = (point) {0, 0};
			ay[i] = (point) {0, 0};
			bx[i] = (point) {0, 0};
			by[i] = (point) {0, 0};
		}
		FFTinit();
		FFT(ax, 1), FFT(ay, 1), FFT(bx, 1), FFT(by, 1);
		static point x[MAXN], y[MAXN], z[MAXN];
		for (int i = 0; i < N; i++) {
			x[i] = ax[i] * bx[i];
			y[i] = ax[i] * by[i] + ay[i] * bx[i];
			z[i] = ay[i] * by[i];
		}
		FFT(x, -1), FFT(y, -1), FFT(z, -1);
		auto num = [&] (point x) {
			return (long long) (x.x + 0.5) % P;
		};
		for (int i = 0; i < N; i++) {
			int res = num(z[i]);
			res = (32768ll * res + num(y[i])) % P;
			res = (32768ll * res + num(x[i])) % P;
			c[i] = res;
		}
	}
}
int n, k, l, x, y, P, g, w;
int mat[3][3], f[MAXN], h[MAXN], res[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;
}
bool PrimitiveRoot(int g) {
	int phi = P - 1;
	for (int i = 2; i * i <= phi; i++)
		if (phi % i == 0) {
			while (phi % i == 0) phi /= i;
			if (power(g, (P - 1) / i) == 1) return false;
		}
	if (phi != 1 && power(g, (P - 1) / phi) == 1) return false;
	return true;
}
int pmat() {
	static int cur[3][3], res[3][3], tmp[3][3];
	for (int i = 0; i <= n - 1; i++)
	for (int j = 0; j <= n - 1; j++) {
		res[i][j] = i == j;
		cur[i][j] = mat[i][j] + (i == j);
		if (cur[i][j] >= P) cur[i][j] -= P;
	}
	int lft = l;
	for (int p = 1; lft != 0; p <<= 1) {
		if (lft & p) {
			lft ^= p;
			for (int i = 0; i <= n - 1; i++)
			for (int j = 0; j <= n - 1; j++) {
				ll tres = 0;
				for (int k = 0; k <= n - 1; k++)
					tres += 1ll * res[i][k] * cur[k][j];
				tmp[i][j] = tres % P;
			}
			memcpy(res, tmp, sizeof(res));
		}
		for (int i = 0; i <= n - 1; i++)
		for (int j = 0; j <= n - 1; j++) {
			ll tres = 0;
			for (int k = 0; k <= n - 1; k++)
				tres += 1ll * cur[i][k] * cur[k][j];
			tmp[i][j] = tres % P;
		}
		memcpy(cur, tmp, sizeof(cur));
	}
	return res[x][y];
}
int main() {
	read(n), read(k), read(l), read(x), read(y), read(P), x--, y--;
	g = 2; while (!PrimitiveRoot(g)) g++; w = power(g, (P - 1) / k);
	for (int i = 0; i <= n - 1; i++)
	for (int j = 0; j <= n - 1; j++)
		read(mat[i][j]);
	for (int i = 0, j = 2 * k - 2; i <= 2 * k - 2; i++, j--) {
		f[i] = power(w, P - 1 - j * (j - 1ll) / 2 % (P - 1));
		if (i <= k - 1) h[i] = 1ll * pmat() * power(w, i * (i - 1ll) / 2 % (P - 1)) % P;
		for (int a = 0; a <= n - 1; a++)
		for (int b = 0; b <= n - 1; b++)
			mat[a][b] = 1ll * mat[a][b] * w % P;
	}
	AnyModuloFFT :: times(f, h, res, P, 2 * k - 2);
	int invk = power(k, P - 2);
	for (int i = 0, j = 2 * k - 2; i <= k - 1; i++, j--)
		writeln(1ll * res[j] * power(w, i * (i - 1ll) / 2 % (P - 1)) % P * invk % P);
	return 0;
}
【Luogu5293】[HNOI2019] 白兔之舞
Yu Gi Oh!
04-09 518
题目链接 题目描述 略 Sol 考场上暴力 O(L)O(L)O(L) 50分真良心。 简单的推一下式子，对于一个 t 来说，答案就是： ∑i=0L[k∣(i−t)](Li)F(i)\sum_{i=0}^{L} [k|(i-t)] {L\choose i}F(i)i=0∑L​[k∣(i−t)](iL​)F(i) 就是对于所有 mod k 的结果是 t 的 i 的后面那一坨东西的和。 F(i)F(i...
【LOJ3057】「HNOI2019」校园旅行
cz_xuyixuan的博客
04-24 547
【题目链接】 点击打开链接 【思路要点】 枚举回文中心向外搜索，可以得到一个 O(N2+M2)O(N^2+M^2)O(N2+M2) 的 BFSBFSBFS 做法。 将边分为两类，一类为同色边，一类为异色边，显然每一步，向外搜索的两个方向 (x,y)(x,y)(x,y) 经过的边应当类型相同。 注意到仅在某一类边上搜索时，我们可以采取如下策略：令 xxx 向目标前进， yyy 在相邻的两...
[HNOI2019]白兔之舞
weixin_33975951的博客
05-07 124
memset0 多合一无聊题 mod k=t，并且k是p-1的约数 单位根反演石锤了。 所以直接设f[i]表示走i步的方案数， 然后C(L,i)分配位置，再A^i进行矩乘得到f[i] 变成生成函数F(x)=∑f[i]=(A*x+I)^L 求指数mod k=t的系数的和 偏移之后，进行单位根反演 对于t都要求？ NTT 然后WA了 因为要任意模数NTT，。。。。 然后...
「HNOI 2019」白兔之舞
weixin_33800463的博客
04-08 121
一道清真的数论题 LOJ #3058 Luogu P5293 题解 考虑$ n=1$的时候怎么做 设$ s$为转移的方案数 设答案多项式为$\sum\limits_{i=0}^L (sx)^i\binom{L}{i}=(sx+1)^L$ 答案相当于这个多项式模$ k$的各项系数的和 发现这和LJJ学二项式定理几乎一模一样 我上一题的题解 然而直接搞是$ k^2$的，无法直接...
Loj 3058. 「HNOI2019」白兔之舞
weixin_34232744的博客
04-12 138
Loj 3058. 「HNOI2019」白兔之舞 题目描述 有一张顶点数为 \((L+1)\times n\) 的有向图。这张图的每个顶点由一个二元组 \((u,v)\) 表示 \((0\le u\le L,1\le v\le n)\)。这张图不是简单图，对于任意两个顶点 \((u_1,v_1),(u_2,v_2)\)，如果 \(u_1<u_2\)，则从 \((u_1,v_1)\) 到 \(...
【LOJ3059】「HNOI2019」序列
cz_xuyixuan的博客
04-24 1183
【题目链接】 点击打开链接 【思路要点】 没有修改的做法在《IOI2018中国国家候选队论文集——浅谈保序回归问题》中有所介绍。 具体做法如下： (1)(1)(1) 、注意到若所有 BiB_iBi​ 均相等，最小化 ∑i=1N(Ai−B)2=∑i=1NAi2−2AiB+B2\sum_{i=1}^{N}(A_i-B)^2=\sum_{i=1}^{N}A_i^2-2A_iB+B^2∑i=1...
【LOJ3073】「2019 集训队互测 Day 2」序列
cz_xuyixuan的博客
05-02 869
【题目链接】 点击打开链接 【思路要点】 不妨令 MMM 为 222 的次幂。 可以发现题目中的初始序列满足：序列的前一半或后一半完全相同，序列的另一半递归地满足该性质。 并且，两个满足该性质的序列卷积后同样满足该性质。 因此，可以用 O(LogM)O(LogM)O(LogM) 个数描述序列，并在 O(LogM)O(LogM)O(LogM) 的时间复杂度内完成卷积。 时间复杂度 O(N...
【LOJ3082】「2019 集训队互测 Day 5」小水题
cz_xuyixuan的博客
05-07 1042
【题目链接】 点击打开链接 【思路要点】 考虑一次钻孔操作对全局水位的影响。 注意到一次钻孔操作只会使得水向某一方向流动，以下认为水向右流动进行讨论。 定义函数 querySum(l,r)querySum(l,r)querySum(l,r) 表示区间 [l,r][l,r][l,r] 的水位和。 定义函数 queryLeft(l,r,h)queryLeft(l,r,h)queryLeft...
【LOJ3180】「IOI2019」天桥
cz_xuyixuan的博客
08-20 1536
【题目链接】 点击打开链接 【思路要点】 类似于子任务 (2)(2)(2) 的解法，考虑找到一些关键点，向四个方向最近的点连边，并在这些点上运行最短路算法。 首先考虑 s=0,t=N−1s=0,t=N-1s=0,t=N−1 的情况，此时初始和目标点均不严格在任意一座天桥的下方，不难证明，此时我们永远不会向左行走。 考虑一座横坐标跨度为 [l,r][l,r][l,r] 的天桥会被如何经过...
【LOJ3058】【HNOI2019】白兔之舞
sadnohappy
04-13 295
Description https://loj.ac/problem/3058 Solution 首先答案长这样子：anst=∑i=0L[k∣(i−t)]Ai(Li)ans_t=\sum\limits_{i=0}^L[k|(i-t)]A^i\binom{L}{i}anst​=i=0∑L​[k∣(i−t)]Ai(iL​)，AAA是读入的矩阵，最后取(x,y)(x,y)(x,y)的值 套上单位根反演就...
HNOI2019 白兔之舞 dance
weixin_30770783的博客
04-09 169
HNOI2019 白兔之舞 dance 显然\(n=3\)就是\(n=1\)的扩展版本，先来看看\(n=1\)怎么做。 令\(W=w[1][1]\)，显然答案是：\(ans_t=\sum_{i\mod k=t}^{L}W^i\binom{L}{i}\) \(=\sum_{i=0}^{L}[k|(i-t)]W^i\binom{L}{i}\) 这时用一个单位根反演。 回顾一下，单位根是fft时用到的...
【LOJ #3058】「HNOI2019」白兔之舞（单位根反演+矩阵快速幂+MTT）
Stargazer的博客
12-28 247
传送门 首先可以发现 走iii步，从xxx到yyy的方案很好算 只需要把方案矩阵AAA求A[x,y]iA^i_{[x,y]}A[x,y]i​再乘一个(Li){L\choose i}(iL​)就可以了 那么对于每个ttt的答案就是 ∑m=0L[k∣m−t](A[x,y]m(Lm))\sum_{m=0}^L[k|m-t](A^m_{[x,y]}{L\choose m})∑m=0L​[k∣m−t](A[...
[HNOI2019]白兔之舞(MTT+单位根反演)
04-23 161
[HNOI2019]白兔之舞 \[ 注:W是一个矩阵,表示题中w[i][j],下列式子中的W做加减法时表示W_{xy}\\ ans_t=\Sigma_{i\ mod\ k=t}^L{L \choose i}W^i\\ =\Sigma_{i=0}^L[k|(i-t)]{L \choose i}W^i\\ =\Sigma_{i=0}^L\frac1k\Sigma_{j=0}^{k-1}(\o...
【loj3058】【hnoi2019】白兔之舞
weixin_34249367的博客
04-28 83
题意 有一个\((L+1)*n\) 的网格图，初始时白兔在\((0,X)\) , 每次可以向横坐标递增，纵坐标随意的位置移动，两个位置之间的路径条数只取决于纵坐标，用\(w(i,j)\) 表示，如果要求白兔停下的点纵坐标为\(Y\) 依次输出移动的步数对\(k\) 取模为 $0 - k -1 $的方案数； \(p\)为质数且$10^8 \lt p \lt 2^{30} , 1 \le n \le ...
loj#3058 「HNOI2019」白兔之舞 MTT+单位根反演+矩阵快速幂
olahiuj的博客
04-11 260
Description 自己看吧 Solution 考场上只想到了n=1的O(L)做法，死于姿势不足不会单位根反演.jpg 正解是真的难写，套一个MTT是最恶心的。。 考虑20分怎么做，一个比较naive的dp就是f[i,j]表示走了i步到了第j列的方案数，nL2转移就可以成功签到了 考虑n=1的情况怎么做。当w=1时，可以发现走恰好d步的方案就是一个组合数，那么答案实际上就是f(d)=∑i=...
LOJ#3058. 「HNOI2019」白兔之舞
04-11 256
LOJ#3058. 「HNOI2019」白兔之舞 https://loj.ac/problem/3058 分析 不妨设转移矩阵为\(f\)，且用\(f^i\)来表示走\(i\)步从\(x\)到\(y\)的方案数。 跳了\(i\)步的方案数为\(\sum\limits_{j=i}^L\binom{j-1}{i-1}=\binom{L}{i}\) 直接推式子\(ans_t=\sum\lim...
loj6279
最新发布
06-21
我们正在查询与LOJ6279相关的编程题目或解决方案。LOJ（LibreOJ）是一个在线的评测系统，主要收录算法竞赛题目。根据题号6279，我们需要确定该题目的具体内容。由于我无法直接访问网络，我将基于已知信息进行推理。...
评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值