几何 - 组合数学 - 分治FFT - 分块

最新推荐文章于 2021-03-27 10:19:26 发布

Mys_C_K

最新推荐文章于 2021-03-27 10:19:26 发布

阅读量436

点赞数

CC 4.0 BY-SA版权

分类专栏：组合计数分治分块快速傅里叶变换

本文链接：https://blog.youkuaiyun.com/Mys_C_K/article/details/82958276

组合计数同时被 3 个专栏收录

35 篇文章

订阅专栏

快速傅里叶变换

23 篇文章

订阅专栏

分块

13 篇文章

订阅专栏

本文探讨了一个关于摧毁正四面体的数学问题，通过组合数学与快速傅里叶变换（FFT）等高级算法解决。核心在于计算特定条件下摧毁正四面体的方案数量，涉及组合数、多项式乘法及模运算。文章提供了详细的算法实现思路与代码示例。

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

题目大意：
有n个正四面体，第k个边长是k个木棍。想要摧毁第k个正四面体，当且进当移除了至少k个正四面体，以及四个结点至少是两条边的段点。现在要摧毁至少m个正四面体，问方案数，不考虑顺序和空间同构， $n≤60000,p=105+3n\le60000,p=10^5+3$ 。p为模数，5组询问。
题解：
首先摧毁第 $k (k > 1)$ 个的方案数是：
$ansk=∑t=04(4t)3t∑i=k−t6k−12(6k−12i)ans_k=\sum_{t=0}^4\binom 4t3^t\sum_{i=k-t}^{6k-12}\binom{6k-12}{i}$
考虑只要知道 $t = 0$ 的情况就能 $O (1)$ 算出其余情况，后办部分是在求组合数后缀和，转为2的幂次减去前缀和，问题转问求：
$bk=S(6k,k+1)∑i=0k+1(6ki)b_k=S(6k,k+1)\sum_{i=0}^{k+1}\binom{6k}{i}$
其实这非常好做，因为： $S(n,k)=2S(n−1,k)−(n−1k)S(n,k)=2S(n-1,k)-\binom{n-1}{k}$ ，因此你可以爆算6次转移。
算出来这个东西之后答案就是 $∑k=mn[xk]∏i=1n(1+ansix)\sum_{k=m}^n\left[x^k\right]\prod_{i=1}^n (1+ans_ix)$ ，这个直接分治FFT即可。
但是由于本人太懒，加上本题时限较大，所以写了个分块优化空间的 $O(n2)O\left(n^2\right)$ 做发并得到了 $100pts100\mathrm{pts}$ 的好成绩。

#include<bits/stdc++.h>
#define rep(i,a,b) for(int i=a;i<=b;i++)
#define lint long long
#define mod 100003
#define N 60060
#define BC N
#define gc getchar()
#define debug(x) cerr<<#x<<"="<<x
#define sp <<" "
#define ln <<endl
using namespace std;
inline int inn()
{
	int x,ch;while((ch=gc)<'0'||ch>'9');
	x=ch^'0';while((ch=gc)>='0'&&ch<='9')
		x=(x<<1)+(x<<3)+(ch^'0');return x;
}
int f[10],ans[N],mi[6*N],fac[6*N],facinv[6*N],b[N],fk[10];
inline int fast_pow(int x,int k,int ans=1)
{
	for(;k;k>>=1,x=(lint)x*x%mod)
		if(k&1) ans=(lint)ans*x%mod;
	return ans;
}
inline int prelude(int n)
{
	int m=min(n,mod-1);
	for(int i=fac[0]=1;i<=m;i++) fac[i]=(lint)fac[i-1]*i%mod;
	facinv[m]=fast_pow(fac[m],mod-2);
	for(int i=m-1;i>=0;i--) facinv[i]=facinv[i+1]*(i+1ll)%mod;
	for(int i=mi[0]=1;i<=n;i++) mi[i]=mi[i-1]*2%mod;return 0;
}
inline int C(int n,int m)
{
	if(n<0||m<0||n<m) return 0;
	if(n>=mod||m>=mod) return (lint)C(n%mod,m%mod)*C(n/mod,m/mod)%mod;
	return (lint)fac[n]*facinv[m]%mod*facinv[n-m]%mod;
}
vector<int> pf[N],s[BC],v;lint w[N];int L[BC],R[BC],bel[N];
inline int tms(vector<int> &a,vector<int> &b,vector<int> &c)
{
	int n=(int)a.size(),m=(int)b.size();
	memset(w,0,sizeof(lint)*(n+m-1));
	for(int i=0;i<n;i++) for(int j=0;j<m;j++)
		w[i+j]+=(lint)a[i]*b[j];c.resize(n+m-1);
	for(int i=0;i<n+m-1;i++) c[i]=w[i]%mod;return 0;
}
int dp[1010][1010];
int main()
{
	f[0]=1,f[1]=12,f[2]=54,f[3]=108,f[4]=81;
	ans[1]=9,ans[2]=243,ans[3]=16224,ans[4]=46489;
	int n=60000;prelude(6*n);b[0]=1;
	rep(k,1,n)
	{
		b[k]=b[k-1];
		rep(i,6*(k-1)+1,6*k) b[k]=b[k]*2-C(i-1,k);
		b[k]=(b[k]%mod+C(6*k,k+1)+mod)%mod;
	}
	rep(k,5,n)
	{
		fk[0]=(mi[6*k-12]-b[k-2]+mod)%mod;
		rep(t,1,4) fk[t]=fk[t-1]+C(6*k-12,k-t);
		rep(t,0,4) ans[k]+=(lint)fk[t]*f[t]%mod;
		ans[k]%=mod;
	}
	int sz=(int)sqrt(n+0.5),bc=(n-1)/sz+1;
	for(int i=1;i<=bc;i++)
	{
		L[i]=R[i-1]+1,R[i]=min(i*sz,n);
		for(int j=L[i];j<=R[i];j++)
		{
			bel[j]=i;
			pf[j].resize(2),pf[j][0]=1,pf[j][1]=ans[j];
			if(j>L[i]) tms(pf[j-1],pf[j],pf[j]);
		}
		s[i]=pf[R[i]];if(i>1) tms(s[i-1],s[i],s[i]);
	}
	int T,x,k,bl;scanf("%d",&T);
	while(T--)
	{
		scanf("%d%d",&x,&k),bl=bel[x];lint ans=0;
		v=pf[x];if(bl>1) tms(s[bl-1],v,v);
		for(int i=k;i<(int)v.size();i++) ans+=v[i];
		printf("%lld\n",ans%mod);
	}
	return 0;
}