poj2888 Magic Bracelet-Burnside

最新推荐文章于 2020-04-17 22:38:30 发布

原创最新推荐文章于 2020-04-17 22:38:30 发布 · 175 阅读

0 ·

CC 4.0 BY-SA版权

ploya 专栏收录该内容

2 篇文章

订阅专栏

本文介绍了一种使用动态规划的矩阵优化方法来解决Polya定理相关问题的创新思路，通过代码实现详细展示了如何利用矩阵乘法和快速幂运算优化计算过程，适用于复杂颜色组合的计数问题。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

这篇博客讲得很好，可以参考一下：https://blog.youkuaiyun.com/maxwei_wzj/article/details/73024349

这是第一次碰到polya定理用动态规划的矩阵优化来做题，前面都是无条件染颜色的基础题。原理引用的博客中讲得非常清楚。

以下是我的草稿，当做是个补充吧

#include<cstdio>
#include<iostream>
#include<cstring>
#define mod 9973
using namespace std;
typedef long long ll;
int m,pr[4000],num,ans;
ll n;
struct matrix
{
	ll s[15][15];
}mp[50];
matrix mult(matrix a,matrix b)
{
	matrix c;
	memset(c.s,0,sizeof(c.s));
	for(int i=1;i<=m;i++)
	{
		for(int j=1;j<=m;j++)
		{
			for(int k=1;k<=m;k++)
			{
				c.s[i][j]=(c.s[i][j]+a.s[i][k]*b.s[k][j])%mod;
			}
		}
	}
	return c;
}
void isprime()
{
	ll k=4000;
	bool prime[4000];
	num=0;
	memset(prime,true,sizeof(prime));
	for(ll i=2;i<k;i++)
	{
		if(prime[i])
		{
			pr[num++]=i;
			for(ll j=i+i;j<k;j+=i)
			{
				prime[j]=false;
			}
		}
	}
}
ll phi(ll x)
{
	ll tp=x;
	for(int i=0;i<num;i++)
	{
		if(x%pr[i]==0)
		{
			tp=tp/pr[i]*(pr[i]-1);
			while(x%pr[i]==0)
			{
				x/=pr[i];
			}
		}
	}
	if(x>1)
		tp=tp/x*(x-1);
	return tp;
}
matrix power(ll x)
{
	int i=0;
  	matrix S;
  	memset(S.s,0,sizeof(S.s));
  	for(int j=1;j<=m;j++) S.s[j][j]=1;
  	while(x)
  	{
  		if(x&1)
  		{
  			S=mult(S,mp[i]);
		}
		x>>=1;
		i++;
	}
	return S;
}
void solve(ll x)
{
  matrix S=power(x);
  ll sum=0;
  for(int i=1;i<=m;i++)
    sum=(sum+S.s[i][i])%mod;
  sum=(sum*phi(n/x))%mod;
  ans=(ans+sum)%mod;
}
void exgcd(ll a,ll b,ll &x,ll &y)
{
  ll x0=1,y0=0,x1=0,y1=1;
  while(b)
  {
    ll tmp,q;
    q=a/b;
    tmp=x0,x0=x1,x1=tmp-q*x1;
    tmp=y0,y0=y1,y1=tmp-q*y1;
    tmp=a,a=b,b=tmp%b;
  }
  x=x0,y=y0;
}
int main()
{
	//FILE *fp=fopen("t.txt","r");
	int cas,k,a,b;
	scanf("%d",&cas);
	isprime();
	while(cas--)
	{
		scanf("%lld%d%d",&n,&m,&k);
		for(int i=1;i<=m;i++)
			for(int j=1;j<=m;j++)
			mp[0].s[i][j]=1;
		while(k--)
		{
			scanf("%d%d",&a,&b);
			mp[0].s[a][b]=mp[0].s[b][a]=0;			
		}
		for(int i=1;i<=35;i++)
			mp[i]=mult(mp[i-1],mp[i-1]);
		ans=0;
		for(ll i=1;i*i<=n;i++)
		{
			if(n%i==0)
			{
				if(i*i==n)
				{
					solve(i);
				}else
				{
					solve(i);
					solve(n/i);
				}
			}
		}
    ll x0,y0;
    exgcd(n,mod,x0,y0);
    x0=(x0%mod+mod)%mod;
    printf("%lld\n",(x0*ans)%mod);
	}
	return 0;
}