「余姚中学 2019 联测 Day 4」随机除法

「余姚中学 2019 联测 Day 4」随机除法

好题，难就难在转移的高位前缀和

首先是一个浅显的\(\text{dp}\)状态，令\(n=\Pi prime_i^{c_i}\)

则状态只跟\(\{c_i\}\)有关，这是一个可重集合，强制定义\(c_i\ge c_{i-1}\)最小表示出所有不同状态

搜索一下\(\text{dp}\)状态，发现只有\(170000\)左右的状态数

直接枚举因数转移复杂度显然是升天的，直接枚举子集状态转移复杂度也很高，并且不好确定系数

所以用一个高位前缀和处理优化枚举因数的转移

再说一遍，是高位前缀和枚举因数的转移，不同的因数可能对应同一个状态

常见的高位前缀和是\(dp_{i,S}=dp_{i-1,S}+dp_{i,S-\{S_i\}}\)

转移具有单调性，对于状态排序之后，定义辅助数组\(dp_{i,j}\)表示对于\(i\)这个状态

它的子集(注意这个子集是未排序的)中和它不同的最低位置\(\ge j\)的总和

计算高位前缀和时，每次转移只会对于一个位置改变

枚举状态时，取得位置是\(j\)，调用时需要排序

而排完序之后\(j\)可能会后移，所以需要定义成\(\ge j\)的，否则会算多

比如转移\((1,1,1)\leftarrow (0,1,1),(1,0,1),(1,1,0)\)

如果定义成\(\le j\)的状态，三种状态转移之后都变成\((1,1,0)\)

原先在这个状态里的三个位置编号是\((0,1,2)\)

如果都去\((1,1,0)\)这个状态里转移过来，原先\((0,1,2)\)对应的下标位置改变，变成

\((1,2,0)\)

\((0,2,1)\)

\((0,1,2)\)

我们访问的时候访问的应该是子状态中不同位置\(\le j\)的总和

而从这个下标改变的状态里转移过来时，原先\(>j\)的下标被移移动进\(\le j\)的范围

再转移就错了

所以正确定义状态之后就可以高位前缀和了

存储和访问状态可以用\(\text{Hash,Trie,int128}\)三种方法存储，\(\text{int128}\)真香啊

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef __int128 Node;
#define rep(i,a,b) for(int i=a,i##end=b;i<=i##end;++i)
#define drep(i,a,b) for(int i=a,i##end=b;i>=i##end;--i)
#define Mod1(x) ((x>=P)&&(x-=P))

const int N=180000,P=1e9+7;

ll qpow(ll x,ll k=P-2) {
	ll res=1;
	for(;k;k>>=1,x=x*x%P) if(k&1) res=res*x%P;
	return res; 
}
Node Max,st[N];
int m,a[100],cnt,dp[N][20],F[N],pri[]={0,2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59,61,67,71,73,79};
char str[30];

int Find(Node s){ return lower_bound(st+1,st+cnt+1,s)-st; }
void Load(Node s){
	rep(i,1,20) a[i]=0;
	for(int i=1;i<=20 && s>1; ++i) while(s%pri[i]==0) s=s/pri[i],a[i]++;
}
void dfs(int p,int lst,Node s) { 
	st[++cnt]=s;
	rep(i,1,lst) {
		if((s*=pri[p])>Max) return;
		dfs(p+1,i,s);
	}
}


int main(){
	freopen("div.in","r",stdin),freopen("div.out","w",stdout);
	Max=1e12,Max=Max*Max;
	dfs(1,100,1);
	fprintf(stderr,"States Number = %d\n",cnt);
	sort(st+1,st+cnt+1);
	rep(i,2,cnt) {
		Node now; Load(now=st[i]);
		int c=1;
		rep(i,1,18) c=1ll*c*(a[i]+1)%P;
		drep(j,18,0) {
			if(!a[j]) dp[i][j]=dp[i][j+1];
			else {
				int k=j;
				while(k>1 && a[k-1]==a[k]) --k;
				int p=Find(now/pri[j]);
				drep(d,j,k) {
					dp[i][d]=dp[i][d+1]+dp[p][d];
					Mod1(dp[i][d]);
				}
				j=k;
			}
		}
		F[i]=(dp[i][1]+c)*qpow(c-1)%P;
		rep(j,1,18) dp[i][j]+=F[i],Mod1(dp[i][j]);
	}
	while(~scanf("%s%d",str,&m)) {
		Node n=0;
		for(int i=0;str[i];++i) n=n*10+(str[i]-'0');
		rep(i,1,m) {
			int x;scanf("%d",&x); a[i]=0;
			while(n%x==0) n=n/x,a[i]++;
		}
		sort(a+1,a+m+1,greater <int> ());
		n=1;
		rep(i,1,m) rep(j,1,a[i]) n=n*pri[i];
		printf("%d\n",F[Find(n)]);
	}
}