[LuoguP5319] [BJOI2019] 奥术神杖 (01分数规划+AC自动机+dp)

[LuoguP5319] [BJOI2019] 奥术神杖 (01分数规划+AC自动机+dp)

题面

神杖上从左到右镶嵌了(n)颗奥术宝石，奥术宝石一共有 1010 种，用数字 0123456789 表示。有些位置的宝石已经残缺，用 . 表示，你需要用完好的奥术宝石填补每一处残缺的部分（每种奥术宝石个数不限，且不能够更换未残缺的宝石）。古老的魔法书上记载了 (m) 种咒语 ((S_i,V_i))，其中 (S_i)是一个非空数字串，(V_i)是这种组合能够激发的神力。

神杖的初始神力值 (mathrm{Magic} = 1)，每当神杖中出现了连续一段宝石与$ S_i$ 相等时，神力值 $mathrm{Magic} $ 就会乘以$ V_i(。但神杖如果包含了太多咒语就不再纯净导致神力降低：设) c$ 为神杖包含的咒语个数（若咒语类别相同但出现位置不同视为多次），神杖最终的神力值为 (sqrt[c]{mathrm{Magic} })。（若 (c=0) 则神杖最终神力值为(1)。）

你需要使修复好的神杖的最终的神力值最大，输出任何一个解即可。

分析

按照套路，把最终答案两边取(ln),每个(V_i)选和不选记为(x_i)

[ln sqrt[c]{Magic}=frac{1}{c} ln prod V_ix_i =frac{1}{c} sum x_i(ln V_i) ]

那么把(V_i)取(ln)之后就变成了一个经典的01分数规划问题。我们可以二分答案(mid),考虑实际值比当前二分值更大的情况,即(frac{1}{c} sum x_i(ln V_i)>mid), 可以转化为选取一些(V_i)来最大化(sum (ln V_i-mid)) ,如果这个最大值大于0,就说明实际值比当前二分值更大。

如何求(max(sum (ln V_i-mid)))?看到多模式串匹配，我们想到在AC自动机上DP. 设(dp_{i,j})表示序列的前(i)位在AC自动机上匹配到节点(j)时的最大值。那么容易写出状态转移方程:

[dp_{i+1,delta(j,k)}=max(dp_{i,j}+cost(delta(j,k))) (k in {0 ,1,2 cdots 9 }) ]

其中(delta(j,k))表示(j)节点通过字符(k)转移后得到的新节点。(cost(x))表示匹配到(x)节点对答案的贡献。如果原串匹配了节点(x)对应的串,答案将增加(V_i-mid).但是,匹配了(x)之后还可能匹配(x)的后缀,因此答案增加的值应该是fail树上(x)到初始节点的路径上所有点的(V_i-mid)之和。每次算一遍(cost)显然会超时,所以我们在BFS预处理时计算(x)到初始节点的路径上所有点的(ln V_i)之和(val_x)和点的个数(sz_x). 这样(cost(delta(j,k))=val_{delta(j,k)}-sz_{delta(j,k)} imes mid)

为了输出答案，还要在dp过程中记录最优转移。

代码

#include<iostream>
#include<cstdio>
#include<cstring>
#include<queue>
#include<cmath>
#include<algorithm>
#define maxn 2000
#define maxc 10
#define INF 1e18
#define eps 1e-9
#define maxln 100
using namespace std;
int n,m;
char str[maxn+5],tmp[maxn+5];
double w[maxn+5];

struct AC_automaton {
	int ch[maxn+5][maxc];
	int sz[maxn+5];//记录x在fail树上到根的链上点的个数
	double val[maxn+5];//记录链上权值和
	int fail[maxn+5];
	int ptr;
	void insert(char *s,double v) {
		int x=0;
		int len=strlen(s+1);
		for(int i=1; i<=len; i++) {
			int c=s[i]-'0';
			if(!ch[x][c]) ch[x][c]=++ptr;
			x=ch[x][c];
		}
		sz[x]++;
		val[x]+=v;
	}
	void get_fail() {
		queue<int>q;
		for(int i=0; i<maxc; i++) if(ch[0][i]) q.push(ch[0][i]);
		while(!q.empty()) {
			int x=q.front();
			q.pop();
			sz[x]+=sz[fail[x]];
			val[x]+=val[fail[x]];
			for(int i=0; i<maxc; i++) {
				if(ch[x][i]) {
					fail[ch[x][i]]=ch[fail[x]][i];
					q.push(ch[x][i]);
				} else ch[x][i]=ch[fail[x]][i];
			}
		}
	}
} T;
double dp[maxn+5][maxn+5];
pair<int,int>res[maxn+5][maxn+5];//记录dp[i][j]由哪个转移过来
void print(int i,int j) {
	if(i==0) return;
	print(i-1,res[i][j].first);
	putchar(res[i][j].second+'0');
}
bool check(double mid,int type) {
	for(int i=0;i<=n;i++){
		for(int j=0;j<=T.ptr;j++){
			dp[i][j]=-INF;
			if(type) res[i][j]=make_pair(0,0);
		}
	}
	dp[0][0]=0;
	for(int i=0; i<n; i++) {
		for(int j=0; j<=T.ptr; j++) {
			if(str[i+1]=='.') { //可以选
				for(int k=0; k<maxc; k++) {
					int nex=T.ch[j][k];
					if(dp[i+1][nex]<dp[i][j]+T.val[nex]-mid*T.sz[nex]) { //01分数规划,每个点变成w[i]-mid
						dp[i+1][nex]=dp[i][j]+T.val[nex]-mid*T.sz[nex];
						if(type) res[i+1][nex]=make_pair(j,k);
					}
				}
			} else {
				int k=str[i+1]-'0';
				int nex=T.ch[j][k];
				if(dp[i+1][nex]<dp[i][j]+T.val[nex]-mid*T.sz[nex]) { //01分数规划,每个点变成w[i]-mid
					dp[i+1][nex]=dp[i][j]+T.val[nex]-mid*T.sz[nex];
					if(type) res[i+1][nex]=make_pair(j,k);
				}
			}
		}
	}
	double ans=0;
	int id=0;
	for(int j=0; j<=T.ptr; j++){
		if(dp[n][j]>ans){
			id=j;
			ans=dp[n][j];
		}
	}
	if(type) print(n,id);
	return ans>=eps;
}

int main() {
	scanf("%d %d",&n,&m);
	scanf("%s",str+1);
	for(int i=1;i<=m;i++){
		scanf("%s",tmp+1);
		scanf("%lf",&w[i]);
		w[i]=log(w[i]);
		T.insert(tmp,w[i]);
	}
	T.get_fail();
	double l=0,r=maxln,mid,ans=0;
	while(fabs(r-l)>eps){
		mid=(l+r)/2;
		if(check(mid,0)){
			ans=mid;
			l=mid;
		}else r=mid;
	}
//	printf("debug: %lf
",exp(ans));
	check(ans,1);
}