BZOJ3557: [Ctsc2014]随机数

BZOJ上数据范围没写清楚,实际数据范围如这个表格所示。


把$M_n$和$x$看成$mathbb{F}_2$下的$m$维的列向量,则有转移矩阵 $$A = egin{pmatrix} 0 & 0 & cdots & 0 & -x_0\ 1 & 0 & cdots & 0 & -x_1\ 0 & 1 & cdots & 0 & -x_2\ vdots & vdots & ddots & vdots & vdots\ 0 & 0 & cdots & 1 & -x_{m-1}\ end{pmatrix}$$ 该矩阵为多项式 $$p(z) = z^m+sum_{i=0}^{m-1}x_iz^i$$ 的友阵。因此,若把$M_n$看成$F$中的元素,则有 $$M_n(z) = z^nM_0(z)$$ 其中$F = mathbb{F}_2[z]/(p(z))$为模$p(z)$意义下的$mathbb{F}_2$上的多项式构成的环。当然,也可以直接通过观察得出这个结论。注意到$n < m$时$F$中的$z^n$与$mathbb{F}_2[z]$中的相同,于是可以以复杂度 $$Oigl(mlog mmax(1,log k-log m)igr)$$ 解决第一种情况。

所有$2^m-1$个非零元素都能被生成,是数列在这些元素中等概率取值的必要条件。该条件成立时有 $$egin{aligned} 2^m-1 &= |mathopen{<}zmathclose{>}M_0(z)| = |mathopen{<}zmathclose{>}|\ &le |F^ imes| le |F setminus {0}| = 2^m-1 end{aligned}$$ 其中$F^ imes$为$F$中所有有乘法逆元的元素构成的乘法群。由此得$mathopen{<}zmathclose{>} = F^ imes = F setminus {0}$,那么$F$为大小为$2^m$的有限域,且$forall f(z) in F^ imes$,有$(f(z))^{2^m} = f(z)$,故 $$egin{aligned} M_k(z) &= (z^{k cdot 2^l})^{2^{m-l}}M_0(z)\ &= Bigl(M_{k cdot 2^l}(z)igl(M_0(z)igr)^{-1}Bigr)^{2^{m-l}}M_0(z)\ &= igl(M_{k cdot 2^l}(z)igr)^{2^{m-l}}Bigl(igl(M_0(z)igr)^{2^{m-l}}Bigr)^{2^l-1}\ end{aligned}$$ 从而可以以复杂度$O(m^2log m)$解决第二种情况。

由上述推导可知,$x$是好的等价于$F$是一个域,且$z$是$F^ imes$的生成元。$F$是一个域当且仅当$p(z)$是$mathbb{F}_2$上的既约多项式。由此,或许可以用某些方法生成第二种情况的数据。


这题挺好写的,就一多项式除法就没了,不知道为啥没几个人写。

#include<algorithm>
#include<vector>
#include<cstdio>
#define RAN(a)a.begin(),a.end()
using namespace std;
typedef unsigned long long u64;
typedef unsigned u32;
namespace num{
	const u32 p=1811939329;
	const u32 g=13;
	inline u32 imod(u32 a){
		return a<p?a:a-p;
	}
	inline u32 ipow(u32 a,u32 n){
		u32 s=1;
		for(;n;n>>=1){
			if(n&1)
				s=(u64)s*a%p;
			a=(u64)a*a%p;
		}
		return s;
	}
	inline u32 inv(u32 a){
		return ipow(a,p-2);
	}
}
using namespace num;
class poly_base{
public:
	vector<u32>a;
	poly_base(){}
	explicit poly_base(int n):a(n){}
	void fft(int n,bool f){
		a.resize(n);
		if(n<=1)
			return;
		for(int i=0,j=0;i<n;++i){
			if(i<j)
				std::swap(a[i],a[j]);
			int k=n>>1;
			while((j^=k)<k)
				k>>=1;
		}
		vector<u32>w(n/2);
		w[0]=1;
		for(int i=1;i<n;i<<=1){
			for(int j=i/2-1;~j;--j)
				w[j<<1]=w[j];
			int m=(p-1)/i/2;
			u64 s=ipow(g,f?p-1-m:m);
			for(int j=1;j<i;j+=2)
				w[j]=s*w[j-1]%p;
			for(int j=0;j<n;j+=i<<1){
				u32*b=&a[0]+j,*c=b+i;
				for(int k=0;k<i;++k){
					u32 v=(u64)w[k]*c[k]%p;
					c[k]=imod(b[k]+p-v);
					b[k]=imod(b[k]+v);
				}
			}
		}
	}
	void dft(int n){
		fft(n,0);
	}
	void idft(){
		int n=a.size();
		fft(n,1);
		u64 f=inv(n);
		for(int i=0;i<n;++i)
			a[i]=f*a[i]%p;
	}
};
class poly:public poly_base{
public:
	poly(){}
	explicit poly(int n):poly_base(n){}
	u32&operator[](int i){
		return a[i];
	}
	const u32&operator[](int i)const{
		return a[i];
	}
	int size()const{
		return a.size();
	}
	void swap(poly&b){
		a.swap(b.a);
	}
	static int len(int n){
		while(n^n&-n)
			n+=n&-n;
		return n;
	}
	void fix(){
		while(size()&&!a.back())
			a.pop_back();
	}
	void mul(poly b){
		fix();
		b.fix();
		int n=len(size()+b.size()-1);
		dft(n);
		b.dft(n);
		for(int i=0;i<n;++i)
			a[i]=(u64)a[i]*b[i]%p;
		idft();
		for(int i=0;i<n;++i)
			a[i]&=1;
		fix();
	}
	void sqr(){
		fix();
		int n=len(2*size()-1);
		dft(n);
		for(int i=0;i<n;++i)
			a[i]=(u64)a[i]*a[i]%p;
		idft();
		for(int i=0;i<n;++i)
			a[i]&=1;
		fix();
	}
	void mod(int n){
		a.resize(n);
	}
	void inv(int n){
		int m=len(n);
		mod(m);
		vector<u32>b(1,1);
		a.swap(b);
		for(int i=2;i<=m;i<<=1){
			poly c(i);
			for(int j=0;j<i;++j)
				c[j]=b[j];
			sqr();
			mul(c);
			mod(i);
		}
		mod(n);
	}
	void div(poly b){
		fix();
		b.fix();
		int n=size()-b.size()+1;
		if(n<=0){
			a.clear();
			return;
		}
		reverse(RAN(a));
		mod(n);
		reverse(RAN(b.a));
		b.inv(n);
		mul(b);
		mod(n);
		reverse(RAN(a));
		fix();
	}
	void mod(poly b){
		fix();
		b.fix();
		int m=b.size();
		if(size()>=m){
			poly c=*this;
			div(b);
			mul(b);
			mod(m-1);
			for(int i=0;i<m-1;++i)
				a[i]^=c[i];
			fix();
		}
	}
};
struct ano{
	operator u64(){
		u64 x=0;
		signed c=getchar();
		while(c<48)
			c=getchar();
		while(c>47)
			x=x*10+c-48,c=getchar();
		return x;
	}
}buf;
poly basis(int n){
	poly a(n+1);
	a[n]=1;
	return a;
}
int main(){
	int n=buf;
	poly p(n+1);
	p[n]=1;
	for(int i=0;i<n;++i)
		p[i]=buf;
	poly s(n);
	for(int i=0;i<n;++i)
		s[i]=buf;
	if(buf==0){
		u64 k=buf;
		if(k<n){
			s.mul(basis(k));
			s.mod(p);
		}else{
			int i=0;
			while((k&(2<<i)-1)<n)
				++i;
			s.mul(basis(k&(1<<i)-1));
			s.mod(p);
			poly a=basis(1<<i);
			k>>=i;
			for(;k;k>>=1){
				if(k&1){
					s.mul(a);
					s.mod(p);
				}
				if(k>1){
					a.sqr();
					a.mod(p);
				}
			}
		}
	}else{
		int l=buf%n;
		poly a(n);
		for(int i=0;i<n;++i)
			a[i]=buf;
		s.swap(a);
		for(int i=0;i<n-l;++i){
			s.sqr();
			s.mod(p);
			a.sqr();
			a.mod(p);
		}
		for(int i=0;i<l;++i){
			s.mul(a);
			s.mod(p);
			if(i<l-1){
				a.sqr();
				a.mod(p);
			}
		}
	}
	s.mod(n);
	for(int i=0;i<n;++i)
		printf("%d",s[i]);
	puts("");
}

Rewrited @ 姫野星奏生誕祭2019

Edited @ 2019-06-28

原文地址:https://www.cnblogs.com/f321dd/p/5556834.html