HDU 5868 Different Circle Permutation

原题链接为：HDU-5868

题意为一群孩子按圈排座位。可以看成是给一个项链染色，有两种颜色0和1，同时相邻处不能都染1，若两种染色方案旋转后相同则视为同一种方案。
从题意可以明显的看出，这是一道考察Burnside引理的数论题。
首先，变换方式只有旋转，可以得出，答案应为( frac{1}{n}sum_{i=1}^{n}f(gcd(i,n)) )。其中f(x)为表示长度为x的一个循环节的染色方案数。

在本题中，考虑最后一个涂1时，则f(x)=f(x-2)，最后一个涂0时，则f(x)=f(x-1)。所以f(x)=f(x-1)+f(x+2)。换言之，斐波那契数列。但是同时必须注意的一点是，f(1)=1。因为在循环节长度为1时，若涂1则项链每个环都是1，不符合题意。所以f(1)=1，f(2)=3，f(3)=4，f(4)=7…本人的做法是将f(0)设为2。同时需要注意的一点是，当n=1是，需要特判。对于斐波那契数列，理所当然的想到要使用矩阵快速幂。

同时，考虑到n的数值非常大，直接使用( frac{1}{n}sum_{i=1}^{n}f(gcd(i,n)) )计算是不可能的。所以对此式进一步简化为(sum _{d mid n} f(gcd(d,n)) * phi(n/d))，其中( phi() )函数，即欧拉函数直接计算即可。

而对于( frac{1}{n} )，因为10^9+7是质数，直接计算逆元即可。

下面是AC的代码：

  1 #include"cstdio"
  2 #include"cstring"
  3 #include"string"
  4 #include"cstdlib"
  5 #include"vector"
  6 #include"set"
  7 #include"map"
  8 #include"cmath"
  9 using namespace std;
 10 typedef long long LL;
 11 
 12 const LL MOD=1000000007;
 13 LL gcd(LL a, LL b)
 14 {
 15     if(b==0) return a;
 16     return gcd(b,a%b);
 17 }
 18 LL powMod(LL x,LL n,LL mod)
 19 {
 20     LL res=1;
 21     while(n>0)
 22     {
 23         if(n&1) res=res*x % mod;
 24         x=x*x % mod;
 25         n>>=1;
 26     }
 27     return res;
 28 }
 29 const LL MAXR=100;
 30 const LL MAXC=100;
 31 struct Matrix
 32 {
 33     LL m[MAXR][MAXC];
 34     LL r,c;
 35     Matrix()
 36     {
 37         r=0,c=0;
 38         memset(m,0,sizeof(m));
 39     }
 40 };
 41 //若不能相乘，返回空矩阵
 42 Matrix operator * (const Matrix &a,const Matrix &b)
 43 {
 44     Matrix ans;
 45     LL r=a.r,c=b.c,x=0;
 46     if(a.c!=b.r) return Matrix();
 47     x=a.c;
 48     ans.r=a.r;
 49     ans.c=b.c;
 50     for (LL i=1; i<=r; i++)
 51         for (LL j=1; j<=c; j++)
 52         {
 53             ans.m[i][j]=0;
 54             for (LL k=1; k<=x; k++)
 55                 ans.m[i][j]=(ans.m[i][j]+a.m[i][k]*b.m[k][j]) % MOD;
 56         }
 57     return ans;
 58 }
 59 //若不能乘方，返回空矩阵
 60 Matrix operator ^ (Matrix &base,LL pow)
 61 {
 62     if(base.r!=base.c) return Matrix();
 63     LL x=base.r;
 64     Matrix ans;
 65     ans.r=ans.c=x;
 66     for (LL i=1; i<=x; i++) ans.m[i][i]=1;
 67     while (pow)
 68     {
 69         if (pow&1) ans=ans*base;
 70         base=base*base;
 71         pow>>=1;
 72     }
 73     return ans;
 74 }
 75 LL eularPhi(LL n)
 76 {
 77     LL res=n;
 78     for(LL i=2;i*i<=n;i++)
 79     {
 80         if(n%i==0)
 81         {
 82             res=res/i*(i-1);
 83             for(;n%i==0;n/=i);
 84         }
 85     }
 86     if(n!=1) res=res/n*(n-1);
 87     return res;
 88 }
 89 //拓展欧几里得算法
 90 //求ax+by=gcd(a,b)的一组解
 91 //其他解为x=x0+kb',y=y0-ka'
 92 //a'=a/gcd(a,b),b'=b/gcd(a,b)
 93 LL extgcd(LL a,LL b,LL &x,LL &y)
 94 {
 95     LL d=a;
 96     if(b!=0)
 97     {
 98         d=extgcd(b,a%b,y,x);
 99         y-=(a/b)*x;
100     }
101     else { x=1; y=0; }
102     return d;
103 }
104 //求逆元
105 //a和ｍ应该互质
106 LL modInverse(LL a,LL m)
107 {
108     LL x,y;
109     extgcd(a,m,x,y);
110     return (m+x%m)%m;
111 }
112 LL f(LL x)
113 {
114     Matrix a,b;
115     a.r=1; a.c=2; a.m[1][1]=2; a.m[1][2]=1;
116     b.r=2; b.c=2; b.m[1][1]=0; b.m[1][2]=b.m[2][1]=b.m[2][2]=1;
117     return (a*(b^x)).m[1][1];
118 }
119 int main()
120 {
121 #ifdef LOCAL
122     freopen("in.txt","r",stdin);
123 #endif
124     LL n;
125     while(scanf("%lld",&n)!=EOF)
126     {
127         if(n==1)
128         {
129             printf("%d
",2);
130             continue;
131         }
132         LL ans=0;
133         for(LL i=1;i<=sqrt(n);i++)
134             if(n%i==0)
135             {
136                 LL d=i;
137                 ans=(ans+f(d)*eularPhi(n/d))%MOD;
138                 d=n/i;
139                 if(d*d==n) continue;
140                 ans=(ans+f(d)*eularPhi(n/d))%MOD;
141             }
142         printf("%lld
",ans*modInverse(n,MOD)%MOD);
143     }
144     return 0;
145 }