【hdu2481】Toy，burnside引理+矩阵乘法

传送门
思路：
快把我做哭了TAT
从昨天上午开始想，搞了一下午有一个点没有想明白
TA爷看过题后想了5min貌似就爆正解了TAT
下面我就来讲一讲~
一开始先想没有置换情况下的方案数
手玩无果后打了表……
然后就是
F(1)=1,F(2)=5,F(n)=3×F(n−1)−F(n−2)+2(n≥3)F(1)=1,F(2)=5,F(n)=3×F(n-1)-F(n-2)+2(n≥3)
mrazer：这不就是“轮状病毒”吗
咦好像真的是啊……
当时那道题我好像也是打表搞得啊……
说说比较好理解的一种证明
设F(n)F(n)是带环的方案数，即n号点可以连到1号点
f(n)f(n)是不带环的方案数，即n号点不能连到1号点
f(n)={12×f(n−1)+∑ni=2f(n−i)(n=0,1)(n≥2)f(n)= \begin{cases} 1&(n=0,1)\\ 2×f(n-1)+\sum^n_{i=2}f(n-i)&(n≥2)\\ \end{cases}
F(n)=∑ni=1i2×f(n−i)F(n)=\sum^n_{i=1}i^2×f(n-i)
先来解释一下F(n)F(n)
考虑把图形外面的n条边想象成分开的一段段
枚举第一个点的段的大小i(1~n)，因为这个段是独立的，想在整体中就必须向中心点连边，所以再枚举这个段向中心点连边的点，数量也是i
剩下的话就是递归问题了，不过要求剩下的段不能成环，就是f(n−i)f(n-i)了
f(n)f(n)的话，考虑前n-1个点已经排好了，加入第n个点使其仍然是树，有三种选择，把n直接连在n-1所在的段上，把n直接与中心点相连，把n连在n-1所在的段上并与中心点相连（也就意味着n-1所在的段形态确定，除n以外不与中心点相连），三种选择的方案分别是f(n−1)f(n-1),f(n−1)f(n-1),∑ni=2f(n−i)\sum^n_{i=2}f(n-i)，最后一种相当于枚举n-1所在段的大小i(1-n)，方案数就取决于剩下的n-1-i个点了
但我们发现上面的式子求解在n≤109n≤10^9下不适用，所以要化简
把f(n) (n≥3)f(n)\ \ (n≥3)化简一下就成了这个样子

f(n)=f(n−1)+∑i=2n−1f(i) ①

f(n)=f(n-1)+\sum^{n-1}_{i=2}f(i)\ \ \ \ \ \ ①

f(n−1)=f(n−2)+∑i=2n−2f(i) ②

f(n-1)=f(n-2)+\sum^{n-2}_{i=2}f(i)\ \ \ \ \ ②
①-②,移项得
f(n)=3f(n−1)−f(n−2) (n≥3)f(n)=3f(n-1)-f(n-2) \ \ (n≥3)
（以下过程我写的有些鬼畜， ~~推荐由结论推出条件~~）

∵∑i=n−2ni2×f(n−i)−3(n−1)2−2(n−2)2−2=0

∵\sum^n_{i=n-2}i^2×f(n-i)-3(n-1)^2-2(n-2)^2-2=0

∴F(n)=∑i=1ni2×f(n−i)−∑i=n−2ni2×f(n−i)+3(n−1)2+2(n−2)2+2

∴F(n)=\sum^n_{i=1}i^2×f(n-i)-\sum^n_{i=n-2}i^2×f(n-i)+3(n-1)^2+2(n-2)^2+2

=∑i=1n−3i2×f(n−i)+3(n−1)2+2(n−2)2+2

=\sum^{n-3}_{i=1}i^2×f(n-i)+3(n-1)^2+2(n-2)^2+2
∵n−i≥3∴f(n−i)=3(n−i−1)−f(n−i−2)∵n-i≥3\\∴f(n-i)=3(n-i-1)-f(n-i-2)

=∑i=1n−3i2×[3f(n−i−1)−f(n−i−2)]+3(n−1)2×f(0)+(n−2)2×[3f(2)−f(1)]+2

=\sum^{n-3}_{i=1}i^2×[3f(n-i-1)-f(n-i-2)]+3(n-1)^2×f(0)+(n-2)^2×[3f(2)-f(1)]+2

=∑i=1n−2i2×[3f(n−i−1)−f(n−i−2)]+3(n−1)2×f(0)+2

=\sum^{n-2}_{i=1}i^2×[3f(n-i-1)-f(n-i-2)]+3(n-1)^2×f(0)+2

=3∑i=1n−1i2×f(n−1−i)−∑i=1n−2i2×f(n−2−i)+2

=3\sum^{n-1}_{i=1}i^2×f(n-1-i)-\sum^{n-2}_{i=1}i^ 2×f(n-2-i)+2

=3F(n−1)−F(n−2)+2

=3F(n-1)-F(n-2)+2
终于搞出了上面的结论了……
对于这种线性递推式就可以矩乘了……
初始矩阵: ⎡⎣⎢F[2]00F[1]00200⎤⎦⎥ \begin{bmatrix} F[2]&F[1]&2\\ 0&0&0\\ 0&0&0\\ \end{bmatrix}
转移矩阵: ⎡⎣⎢3−11100001⎤⎦⎥ \begin{bmatrix} 3&1&0\\ -1&0&0\\ 1&0&1\\ \end{bmatrix}
也有是用行列式、基尔霍夫矩阵直接秒的……蒟蒻不太会，准备找个时间去学习一下O__O
然后就是套burnside引理了
我们要考虑的是不动点的数量，所以前 gcd(i,n)gcd(i,n)个点的连边方式确定下来，整个图的连边方式就确定了
我一开始想的是 f[gcd(i,n)]f[gcd(i,n)]，但这是不对的
因为这些点也是可以看成分段的
最后一个点的下一个点可以看成是第一个点(即最后一段的下一段可以看成是第一段)，实际上这是一个环，所以这 gcd(i,n)gcd(i,n)个点的连接方式与 F[gcd(i,n)]F[gcd(i,n)]是等价的，即合法性等效
所以答案是 ∑ni=1F[gcd(i,n)]\sum^n_{i=1}F[gcd(i,n)]
化简以后就是 ∑ni=1[i|n]F(i)φ(ni)\sum^n_{i=1}[i|n]F(i)\varphi(\frac n i)
万里长征就差最后一步了
除个n就可以了
啥？你说n,m可能不互质？
难道还要补个中国剩余定理……
有一个好用的式子
a/bmodp=(amodbp)/ba/b\mod p=(a\mod bp)/b,其中 b|ab|a
证明也十分简单
设 a=bka=bk
a/bmodp=(bk)/bmodp=kmodpa/b\mod p=(bk)/b\mod p=k\mod p
(amodbp)/b=(bkmodbp)/b=b(kmodp)/b=kmodp(a\mod bp)/b=(bk\mod bp)/b=b(k\mod p)/b=k\mod p
所以运算时模 nmnm，最后除以m就可以了
nm≤1018nm≤10^{18}所以做乘法的时候还要用快速乘……

#include<cstdio>
#include<cmath>
#include<cstring>
#include<iostream>
#define LL long long
using namespace std;
int n,m;
int prime[32000];
bool vis[32000];
LL mo,ans;
LL mul(LL x,LL y)
{LL t=0;bool flag=0;if (y<0) y=-y,flag=1;for (;y;y>>=1,x=(x+x)%mo)if (y&1) t=(t+x)%mo;return flag?-t:t;
}
struct Matrix
{LL a[4][4];void clr(){memset(a,0,sizeof(a));}
}A,B;
Matrix operator *(Matrix A,Matrix B)
{Matrix C;C.clr();for (int i=1;i<=3;++i)for (int j=1;j<=3;++j)for (int k=1;k<=3;++k)C.a[i][j]=(C.a[i][j]+mul(A.a[i][k],B.a[k][j]))%mo;return C;
}
LL cal(int x)
{if (x==1) return 1;if (x==2) return 5;A.clr();B.clr();A.a[1][1]=5;A.a[1][2]=1;A.a[1][3]=2;B.a[1][1]=3;B.a[2][1]=-1;B.a[1][2]=1;B.a[3][1]=1;B.a[3][3]=1;x-=2;for (;x;x>>=1,B=B*B)if (x&1) A=A*B;return A.a[1][1];
}
int phi(int x)
{int ans=x;for (int i=1;prime[i]*prime[i]<=x;++i)if (x%prime[i]==0){ans=ans/prime[i]*(prime[i]-1);while (x%prime[i]==0) x/=prime[i];}if (x>1) ans=ans/x*(x-1);return ans%mo;
}
main()
{int lim=31625;for (int i=2;i<=lim;++i){if (!vis[i])prime[++prime[0]]=i;for (int j=1;j<=prime[0];++j){if (i*prime[j]>lim) break;vis[i*prime[j]]=1;if (i%prime[j]==0) break;}}while (~scanf("%d%d",&n,&m)){mo=1LL*n*m;ans=0;for (int i=1;i*i<=n;++i)if (n%i==0){ans=(ans+mul(phi(n/i),cal(i)))%mo;if (i*i!=n) ans=(ans+mul(phi(i),cal(n/i)))%mo;}printf("%d\n",(ans/n+m)%m);}
}

【hdu2481】Toy，burnside引理+矩阵乘法相关推荐

HDU 2865 Birthday Toy [Polya 矩阵乘法]
传送门题意: 相邻珠子不能相同,旋转等价.$n$个珠子$k$中颜色,求方案数首先中间珠子$k$种选择,$k--$ 如果没有相邻不同的限制,就和$POJ\ 2154$一样了 $|C(f)|=k^{\ ...
POJ 2888 Magic Bracelet ——Burnside引理
[题目分析] 同样是Burnside引理.但是有几种颜色是不能放在一起的. 所以DP就好了. 然后T掉所以矩阵乘法就好了. 然后T掉所以取模取的少一些,矩阵乘法里的取模尤其要注意,就可以了. A掉 ...
快速矩阵乘法的研究——下
快速矩阵乘法的研究本文我们只看一个算法,就是 Coppersmith 和 Winograd 提出的 O ( 2.37 ) O(2.37) O(2.37)的矩阵乘法.请确保前两篇内容已经掌握理解. h ...
Polya定理与Burnside引理及其应用
Polya定理及其应用群与置换群轮换 PolyaPolyaPolya定理 BurnsideBurnsideBurnside引理 PolyaPolyaPolya定理的简单应用 PolyaPolyaP ...
如何在CPU上优化GEMM矩阵乘法
如何在CPU上优化GEMM矩阵乘法 How to optimize GEMM on CPU (TL;DR) TVM 提供抽象接口,允许用户分别描述算法和算法的实现组织(所谓的调度).通常,在高性能调度 ...
CPU的自动调度矩阵乘法
CPU的自动调度矩阵乘法这是一个有关如何对CPU使用自动调度程序的文档. 与依靠手动模板定义搜索空间的基于模板的autotvm不同,自动调度程序不需要任何模板.用户只需要编写计算声明,而无需任何调度 ...
十个利用矩阵乘法解决的经典题目
出自matrix67.com 好像目前还没有这方面题目的总结.这几天连续看到四个问这类题目的人,今天在这里简单写一下.这里我们不介绍其它有关矩阵的知识,只介绍矩阵乘法和相关性质. 不要以为数学 ...
[学习笔记]矩阵乘法及其优化dp
1.定义: $c[i][j]=\sum a[i][k]\times b[k][j]$ 所以矩阵乘法有条件,(n*m)*(m*p)=n*p 即第一个矩阵的列数等于第二个矩阵的行数,否则没有意义. 2.结 ...
ICML 2021：矩阵乘法无需相乘，速度提升100倍，MIT开源最新近似算法
点击上方"视学算法",选择加"星标"或"置顶" 重磅干货,第一时间送达萧箫发自凹非寺量子位报道 | 公众号 QbitAI 在不做乘 ...

【hdu2481】Toy，burnside引理+矩阵乘法

【hdu2481】Toy，burnside引理+矩阵乘法相关推荐

最新文章

热门文章