Invoker

2
3 4
1 2

Case #1: 21
Case #2: 1

Hint
For Case #1: we assume a,b,c are the 3 kinds of elements.
Here are the 21 different arrangements to invoke the skills
/ aaaa / aaab / aaac / aabb / aabc / aacc / abab / 
/ abac / abbb / abbc / abcb / abcc / acac / acbc /
/ accc / bbbb / bbbc / bbcc / bcbc / bccc / cccc /
 

Polya计数。题目可转化为用c种颜色给n个珠子的项链染色，问一共同拥有多少种颜色方案。本题要对结果取模1000000007

1.旋转。

将环顺时针旋转i格后，循环节个数为gcd(n,i)， 染色方案为  ∑c^gcd(n,i)    当中 i=1,2,3,4,....n

2.翻转。

这里也得考虑两种情况。

        当n为奇数时。共同拥有n个循环节个数为（n/2+1)的循环群，还有的资料上说是环的个数为（n/2+1) ，注意这是计算机上的表示，n/2整型相除计算机得到的是整数。事实上应该写成（n+1)/2。。染色方案为  n*c^(n/2+1)

为什么n个循环节个数为（n/2+1)的循环群呢？我的理解是这种，也许不太对。。。

       拿正三角形为例，给它三个顶点染色， 对称轴是一个顶点与其对边终点连线所在的直线，这种直线有3（n=3。即n个顶点) 条，共同拥有3（n)个循环群。

如果第一个顶点在对称轴上，那么第二个顶点经过对称轴翻转肯定和第三个顶点重合，那么 （2,3）是一个循环节。（1）自己是一个循环节，循环节个数为2，即（n+1/2)。

        当n为偶数时，共同拥有n个循环群。当中有n/2个的循环节个数为（n/2 +1)。 有n/2个的循环节个数为（n/2)。

拿正方形为例。四个顶点从左上角顺时针编号1,2,3,4.  

当以1,3顶点连线所在直线为对称轴时（对角的两个顶点）,这样对称轴有2个（n/2)。经过翻转。2,4 重合，1和1重合，3和3重合。那么循环节的个数为3（2,4） （1）（3）, 即(n/2+1)。 染色方案为  (n/2)*c^(n/2+1)

当以两条相对平行的边的中点连线所在直线为对称轴时，比方以线段1,2的中点和3,4的中点连线的所在直线为对称轴。这种对称轴有两个（n/2)，经过翻转。1,2重合，3,4重合，循环节的个数为2。（1,2）（3,4），即（n/2)。，也就是谁和谁重合。谁就和谁在一个循环节里。

染色方案为(n/2)*c^(n/2)

最后累加方案得到ans， 再除以置换群的个数2*n，即 ans/(2*n)%mod即为最后答案。但这里要特别注意，ans是在计算过程中不断取模得到的数。ans,2*n都在模剩余系中，不能直接參与除法计算，由于有公式a*b%mod=(a%mod*b%mod)%mod，除法对取余不满足结合律，a/b!=((a%mod)/(b%mod))%mod ，在计算 ans/(2*n）%mod时。能够转化为 ans*inv(2*n)%mod ,当中  inv(2*n)是2*n关于mod的逆元，保证乘以inv(2*n)和除以 2*n 对于最后的答案取余mod是一样。 

所以如今的问题是如何求一个数关于模P的逆元。

方法1：扩展欧几里得。  ax=1(mod P), gcd(a,p)=1。 当中x为a的逆元，就是我们所求，ax=PY+1,   ax-Py=1,  所以用扩展欧几里得能够求出x。

方法2：费马小定理：  假设模P是素数的话，那么inv(a)=pow(a,p-2)%p; 等式右边用高速幂运算能够得出。

參考：http://www.xuebuyuan.com/1394391.html

代码：

#include <iostream>
using namespace std;
typedef long long LL;
const LL mod=1000000007;
LL c,n;

LL gcd(LL a,LL b)
{
    return b==0?
a:gcd(b,a%b);
}

LL power(LL p,LL n)//高速幂运算
{
    LL ans=1;
    while(n)
    {
        if(n&1)
            ans=ans*p%mod;
        p=p*p%mod;
        n/=2;
    }
    return ans;
}
LL  exgcd(LL a,LL b,LL &x,LL &y)//扩展欧几里得算法。返回a,b的最大公约数,ax+by=gcd(a,b)，x,y为方程的一组解
{
    if(b==0)
    {
        x=1;
        y=0;
        return a;
    }
    long long d=exgcd(b,a%b,x,y);
    long long t=x;
    x=y;
    y=t-a/b*y;
    return d;
}

int main()
{
    int t;cin>>t;
    int cas=1;
    while(t--)
    {
        cin>>c>>n;
        int ans=0;
        for(LL i=1;i<=n;i++)
        {
            ans+=power(c,gcd(n,i));
            ans%=mod;
        }
        if(n&1)
            ans+=(n*power(c,n/2+1))%mod;
        else
            ans+=((n/2*power(c,n/2+1))%mod+(n/2*power(c,n/2))%mod)%mod; //注意mod的位置
        ans%=mod;
        LL x,y;
        exgcd(2*n,mod,x,y);
        //x=power(2*n,mod-2)%mod;//另外一种方法
        x=(x+mod)%mod;
        cout<<"Case #"<<cas++<<": "<<ans*x%mod<<endl;
    }
    return 0;
}