[多校联考2021] 妹子公司

一、题目

题目背景

哥哥 \(\tt zxy\) 的公司有无限个妹子，她们要开演唱会 \(...\)

题目描述

有一个 \(m\) 个数位的数，每个数位的数字 \(\in\{1,2,3,5,7\}\)，最后的数字不能被 \(2,3,5,7,11,47\) 中的任意一个整除，求合法数字个数，答案对质数 \(9973\) 取模。

数据范围

\(m\leq10^9\)

二、解法

首先有一个小优化，由于是否被 \(2,5\) 整除是由末尾数字决定的，所以可以先不考虑这两个除数，单独考虑最后一个数就行了，那么设 \(n=3\cdot 7\cdot 11\cdot 47\)，我们算出模 \(n\) 为 \(x\) 的数字个数就可以得到答案。

然后是套路地数位 \(dp\)，设 \(dp[i][j]\) 为考虑了 \(i\) 个数位，数字模 \(n\) 为 \(j\) 的方案数，每次考虑加入一个数字，你发现第一维太大所以可以用矩阵加速，时间复杂度 \(O(n^3\log m)\)

但是你发现还是会炸成傻逼，谁告诉你加速递推只能用矩阵了？还可以用多项式倍增加速递推：

\[dp[i][(j\cdot 10^{i/2}+k]=\sum dp[\frac{i}{2}][j]\times dp[\frac{i}{2}][k] \]

这个东西不能再像多项式乘法了，只需要设 \(h[i]=\sum_{j\cdot 10^{i/2}=i\bmod n}dp[\frac{i}{2}][j]\)，然后做循环卷积即可，时间复杂度 \(O(n\log n\log T)\)，\(\tt FFT\) 要打递推版不然会凉，写多项式题要注意传进去的数组最好不要随便改。

#include <cstdio>
#include <iostream>
#include <cmath>
using namespace std;
#define db double
const int M = 100005;
const int MOD = 9973;
const db pi = acos(-1.0);
#define ll long long
int read()
{
	int x=0,f=1;char c;
	while((c=getchar())<'0' || c>'9') {if(c=='-') f=-1;}
	while(c>='0' && c<='9') {x=(x<<3)+(x<<1)+(c^48);c=getchar();}
	return x*f;
}
int n,m,k,ans,a[M],r[M],f[M],rev[M];
namespace poly
{
	struct comp
	{
		db x,y;
		comp(db X=0,db Y=0) : x(X) , y(Y) {}
		comp operator + (const comp &R) const {return comp(x+R.x,y+R.y);}
		comp operator - (const comp &R) const {return comp(x-R.x,y-R.y);}
		comp operator * (const comp &R) const {return comp(x*R.x-y*R.y,x*R.y+y*R.x);}
	}A[M],B[M];
	void FFT(int len,comp *a,int fl)
	{
		for(int i=0;i<len;i++)
		{
			rev[i]=(rev[i>>1]>>1)|((len/2)*(i&1));
			if(i<rev[i]) swap(a[i],a[rev[i]]);
		}
		for(int s=2;s<=len;s<<=1)
		{
			int t=s/2;comp w=comp(cos(2*pi/s),sin(2*pi/s)*fl);
			for(int i=0;i<len;i+=s)
			{
				comp x=comp(1,0);
				for(int j=0;j<t;j++,x=x*w)
				{
					comp fe=a[i+j],fo=a[i+j+t];
					a[i+j]=fe+x*fo;
					a[i+j+t]=fe-x*fo;
				}
			}
		}
	}
	void work(int n,int *a,int *b)
	{
		int len=1;while(len<2*n) len<<=1;
		for(int i=0;i<len;i++) A[i]=B[i]=comp(0,0);
		for(int i=0;i<n;i++) A[i]=comp(a[i],0),B[i]=comp(b[i],0);
		FFT(len,A,1);FFT(len,B,1);
		for(int i=0;i<len;i++) A[i]=A[i]*B[i];
		FFT(len,A,-1);
		for(int i=0;i<n;i++) b[i]=0;
		for(int i=0;i<len;i++)
			b[i%n]=(b[i%n]+(ll)((A[i].x/len)+0.5))%MOD;
	}
}
int qkpow(int a,int b)
{
	int r=1;
	while(b>0)
	{
		if(b&1) r=1ll*r*a%n;
		a=1ll*a*a%n;
		b>>=1;
	}
	return r;
}
void fuck(int *a,int *b)
{
	int t=qkpow(10,k);
	for(int i=0;i<n;i++) f[i]=0;
	for(int i=0;i<n;i++) f[1ll*i*t%n]=(f[1ll*i*t%n]+b[i])%MOD;
	poly::work(n,a,f);
	for(int i=0;i<n;i++) b[i]=f[i];
}
signed main()
{
	freopen("c.in","r",stdin);
	freopen("c.out","w",stdout);
	n=3*7*11*47;m=read();m--;k=1;
	r[0]=a[1]=a[2]=a[3]=a[5]=a[7]=1;
	while(m>0)
	{
		if(m&1) fuck(a,r);
		fuck(a,a);
		k<<=1;
		m>>=1;
	}
	int w[10]={1,3,7,11,47};
	for(int i=0;i<n;i++)
		for(int j=0;j<=2;j++)
		{
			int t=(i*10+w[j])%n,f=1;
			for(int k=1;k<=4;k++)
				if(t%w[k]==0) f=0;
			if(f) ans=(ans+r[i])%MOD;
		}
	printf("%d\n",ans);
}